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Mle  Theory 

Discrete  Ordinates  Methods 
Light  Scattering 
Legendre  Polynomials 


ICMtffMM  on  r«««r««  o9^  it  rnmommomw  mtS  tf  M««A 

A study  was  performed  to  determine  the  applicability  of  discrete  ordi- 
nates codes  to  the  transport  of  light  through  particulates.  Light  Is  scat- 
tered with  high  anisotropy  by  particulates.  This  high  anisotropy  requires 
, expansion  In  high  Legendre  order  for  accurate  expression  by  discrete  ordi- 
nates codes. 

Procedures  for  convenient  and  economical  production  of  light  scattering 
data  (Mle  theory)  In  discrete  ordinates  format  were  developed  and^ested.  — ■ 
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One-dlaenslonal  tlae-dependenf  dl89rete  ordinates  calculations  were  performed 
using  the  TDA  (Time  Dependent  ANISN)  code  with  a simple  cloud  model.  Com- 
puter storage  and  time  limitations  are  outlined. 

Extension  to  largest  particle  sizes  and  a second  space  dimension  can 
best  be  accomplished  using  asymmetric  discrete  angle  sets  to  conserve  com- 
puter resources.  Further  economies  might  be  accomplished  by  normalization 
of  low-order  data  for  problems  dominated  by  multiple  scattering.  Relaxation 
of  the  need  to  express  data  In  Legendre  series  could  offer  further  economy, 
but  may  require  significant  development. 
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I.  INTRODUCTION 


Discrete  ordinates  and  other  codes  for  radiative  transfer  calculations 
require  phase  function  Information  to  be  expressed  In  terms  of  Legendre 
polynomials.  A phase  function,  evaluated  by  Mle  theory,  can  be  expressed 
directly  In  Legendre  polynomials  by 

N 

p(cos0)  - I P^(cos0)  (1) 

—ft  li  n ' ' 

n-0 

where  P^(cos0)  Is  the  Legendre  polynomial  of  order  n and  N*  Is  the  maximum 

order  necessary  to  sufficiently  express  the  phase  function  and  a Is  the 

n 

associated  Legendre  coefficient. 

The  phase  function  p(cos0)  Is  a function  of  scattering  particle 
size,  radiation  wavelength,  and  refractive  Index.  Thus,  sets  of  a can 

R 

be  determined  for  each  combination  of  size,  wavelength,  and  refractive 
Index.  Normally,  one  would  wish  to  perform  transport  calculations  for  a 
single  wavelength  and  refractive  Index,  but  over  a range  of  scattering 
particle  sizes.  An  additional  parameter  Is  defined: 

a ■ 2Tr  r/X  (2) 


where 


Cl  ’•  size  parameter 
r - particle  radius 
A ■ radiation  wavelength. 

Then  for  a given  refractive  Index  and  wavelength,  a range  of  particle  sizes 

Is  translated  Into  a range  of  ct-values.  Sets  of  Legendre  coefficients,  a , 

n 

can  then  be  determined  for  each  combination  of  refractive  Index  and  size 
parameter  a. 

a 

Mathematics  tradition  generally  dictates  the  use  of  n«0  for  the  lowest  order 
considered,  but  computer  Indexing  usually  does  not  allow  the  use  of  zero,  so 
equations  related  to  computer  progranmlng  generally  Indicate  n-l,N.  Thus, 
the  N of  a computer-related  equation  may  be  1 greater  than  N In  an  equivalent 
equation  In  traditional  form. 

1 

ss 


Bi 


There  will  generally  be  a different  number  of  particles  for  each 
3 

size,  r.  The  number  per  cm  for  any  r Is  the  particle  size  distribution 
n(r)  such  that  /n(r}dr  Is  the  density  of  particles,  partlcles/cm'^.  In  the 
size  range  Ar .^^Appropriately  averaged  phase  functions  can  be  determined 
by  Integration  of  single  size  phase  functions  over  an  appropriate  distribution 
n(r). 

Transport  codes  require  data  In  the  form  of  Interaction  probability 
per  unit  path  length  (macroscopic  cross  section)  rather  than  the  relative 
fractional  scatter  per  particle  given  by  the  phase  function.  This  Interaction 
probability  can  be  expressed  In  the  same  expansion  form  as  the  phase  function: 

N 

Z (co80)“  y A P (cos6)  (3) 

8 r\  n n 

n*0 


where  the  coefficients, 
®n 


A are  related  to  the  phase  function  coefficients, 
n 


emd  where 


Og(a,B)  - 

o ■ 

e - 


m " 


1 - 


n(r) 


10-8- 


a (a,3)a  (a,e)n(r)dr 
J s n 

Ar 

microscopic  scattering  cross  section 

size  parameter,  2Trr/X 

a(m-lk) 

Index  of  refraction  relative  to  the 
surrounding  medium 

extinction  coefficient  of  the  particle 
material 

Imaginary  coefficient  for  the 
complex  number  m 

particle  number  density  function  for 
particles  of  size  r 
cm^ /micron^ . 


(4) 


A set  of  procedures  for  calculating  Interaction  probability  coeffi- 
cients has  been  produced  and  tested  with  the  Time-Dependent  AMISM  (TDA)  one- 


space-dimension  discrete  ordinates  code  (Ref.  1).  A description  of  the 
procedures  and  recommendations  based  on  the  tests  follows.  The  appendices 
give  utilization  instructions  for  the  computer  codes  used  in  this  study 
and  a file-name  directory  for  the  files  stored  on  the  Xerox  Slgma-9  computer 


II.  APPLICATION  TO  DISCRETE-ORDINATE 

transport  codes 


The  use  of  Legendre  expansions  of  cross  section  Is  widely  practiced 

In  computer  calculations  of  radiation  transport.  Computer  codes  for  large 

radiation  transport  calculations  can  generally  be  grouped  Into  two  categories: 

(1)  Monte  Carlo,  based  on  random  representation  of  natural  processes,  and  (2) 

Boltzmann  balances,  based  on  the  Boltzmann  equation.  Generally,  the  Monte 

Carlo  method  allows  a wide  choice  of  methods  by  which  the  cross-section 

data  can  be  represented.  The  Boltzmann  balance  method  codes,  however, 

usually  use  cross  sections  represented  in  Legendre  pol3momlals.  Prominent 

among  the  Boltzmann  balance  codes  are  those  using  the  S -discrete  ordinate 

n 

numerical  techniques. 

The  S -discrete  ordinate  codes  are  limited  to  either  one  or  two 
n 

spatial  dimensions,  and  the  choice  between  them  depends  on  the  geometry 
of  the  case  being  modeled. 

2.1  Geometry  Limitations 

2.1.1  One  Dimension 

One-dlmenslonal  codes  (ANISN  (Ref.  2),  TDA,  DRF  (see  Ref.  2)) 
are  severely  limited  In  the  extent  to  which  they  can  be  used  to  model 
real  situations  (see  Fig.  1). 

Slab  Geometry:  If  the  scattering  medium  Is  considered  as  an 
Infinite  slab,  with  the  allowed  dimension  the  slab  thickness,  then  the 
source  must  be  an  Infinite  plane  source.  Some  code  versions  allow  a choice 
of  source  directions,  but  all  require  spatial  Infinity. 

Spherical  Geometry:  If  the  scattering  medium  Is  a sphere,  with  the 
allowed  dimension  the  radial  thickness  of  the  Increment,  a point  source  Is 
possible.  A point  source,  however,  must  be  located  at  the  center  of  the 
sphere.  Location  at  any  radius  greater  than  zero  causes  It  to  become  a 
shell  source. 

Cylindrical  Geometry:  If  the  scattering  medium  Is  an  Infinitely  long 
cylinder,  with  the  allowed  dimension  the  radial  thickness  of  the  medium,  the 
source  can  be  an  Infinitely  long  line  source  at  the  center  or  a cylindrical 


shell  source  In  any  radial  Increment.  The  requirement  of  a point  source 

necessarily  limits  one-dlmenslonal  simulation  to  spherical  geometry.  This, 

then,  requires  that  the  scattering  medium  have  curvature  equal  to  the 

distance  from  the  source.  This  further  requires  that  all  of  the  source 

radiation  be  emitted  normal  to  the  Incident  cloud  surface  and  that  the 

detector  position  be  determined  by  the  radial  distance  from  the  source 

only.  The  real  conditions  of  varying  slant  thickness  and  cloud- to-detector 

distance  might  be  approximated  by  choosing  appropriate  combinations  of 

Interval  thickness  and  leakage  angle  (see  Fig.  2).  ' 

This  approximation  is  probably  best  for  relatively  thin  clouds, 
since  there  would  be  fewer  contributions  from  multi-scatters  that  varied 
greatly  from  a path  along  the  common  thickness.  The  approximation  Is 
aided  by  the  general  forward-scattering  character  of  light.  Where  the 
radiation  path  varies  greatly  from  a straight  line  between  the  source  and 
detector,  the  approximate  leakage  angle  6*  Is  much  different  from  the  real 
leakage  angle  6.  If  the  difference  between  them  Is  great,  especially  when 
one  of  them  approaches  0°  (a  90°  scatter  change) , the  approximation  may  be 
poor. 

2.1.2  Two  Dimensions 

The  limitations  described  for  one-dlmenslonal  representations  are 
largely  avoided  by  using  a two-dimensional  representation  (DOT  (Ref.  3), 

TRANZIT  (Ref.  4),  and  TWOTRAN  (Ref.  5).  The  sample  case  (Fig.  1),  can 
be  well-represented  in  cylindrical  geometry  (p,2)  by  describing  a region 
of  cloud  around  a normal  from  the  source  (see  Fig.  3).  The  source  Itself 
can  be  represented  by  Intervals  in  p and  z around  p-0.  Then  the  detector 
response  can  be  determined  by  Integrating  over  all  azimuthal  positions  i 

in  the  upper  surface  r-lntervals.  Intervals  can  be  defined  to  values  of 
r as  large  as  necessary  to  represent  the  detector  response.  The  azimuthal 
symmetry  assumed  in  the  two-dimensional  case  Is  not  a severe  limitation 
as  such  symmetry  can  be  reasonably  assumed. 
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g.  2.  One-dlmenslonal  Representation  of  Sample  Case 


Many  analytical  and  most  numerical  applications  of  the  Boltzmann 


transport  equation  employ  Legendre  polynomials  In  some  form.  The  reason 
for  this  lies  In  the  form  of  the  scattering  term,  which  gives  the  density  of 
particles  (or  photons)  scattered  Into  the  direction  n from  all  other 
directions  U': 


dSl' 

4ir 

where  is  the  flux  in  direction  H'.  The  scattering  cross  section, 

£ Is  then  expanded  In  Legendre  polynomials;  for  Instance,  by: 


(5) 


^ ! (^>  I, 

n“0 


(6) 


where  £ Is  the  nrt-lst  Legendre  coefficient  and  P (p  ) Is  the  n^-lst  Legendre 
n no 

polynomial.  The  argument,  y^.  Is  the  cosine  of  the  angle  between  the  two 
vector  directions,  (i'  and  fl.  This  expansion  Is  very  convenient  for  neutrons 
because.  In  many  cases.  Isotropic  scatter  predominates.  In  Isotropic  scatter, 
the  Initial  and  final  directions  do  not  correlate  and  the  probability  of 
scatter  into  any  unit  solid  angle  Is  a constant.  In  this  case,  only  the 
first  term  of  the  expansion  needs  to  be  used,  where  Is  s constant 

(P^(y^)  ■ 1).  While  this  simplified  case  Is  seldom  exact.  If  often  performs 
well  as  a first  approximation  and  good  results  are  often  obtained  for 
relatively  low  orders  of  expansion  (P^  or  P^) . 

Scattering  of  gamma  rays  Is  not  well-approximated  by  Isotropic 
scattering,  but  Is  of  sufficiently  low  order  that  good  results  are  usually 
obtained  by  calculations  employing  only  a few  more  terms  than  are  used  for 
neutrons.  Light  scattering,  however,  can  be  of  a very  high  order  when 
cross  sections  are  expressed  In  Legendre  polynomials.  Orders  of  several 
hundred  may  be  required  to  express  Mle  theory  phase  functions  for 
relatively  large  aerosols  (a  few  tens  of  microns  In  radius) . 


The  convenience  of  the  Legendre  polynomial  expansion  In  the  scattering 
term  lies  in  the  form  of  the  polynomial  P^(y^).  It  can  be  written  in  terms  of 
the  angles  defining  the  Initial  and  final  dlrectlonSf  and 


where 


and 


P (y  ) = P (H'.n) 
no  n 

» P„(y)P  (y')  + 
n n 


pf(y)  = (-i)^(i-y")®^^ 


2 I fclH  p?(y’)co8ea-i;') 

d®p^(y) 


(7) 


y'  and  y are  the  cosines  of  the  polar  angles  describing  the 
directions  ?7'  and  Q, 

and  ^ are  the  azimuthal  angles  describing  the  directions 
n'  and  n. 


Then  the  scattering  term  can  be  written  in  two  parts: 


1_  J .2ntl.j. 

2Tr  2 ^ n 

n*0 


2-n  1 

r 


dC' 

0 -1 


dy'  P (y)  P (y')  (|>(n*) 
n n 


i?  (IStlxj; 

IT  ^ 2 ^ n 

n-0 


2 IT  1 
f 


d^' 

0 -1 


J 


"(53)1  pJ(y')co8B(c-i:')^(n') 


3-1 


(8) 


(9) 


Assuming  azimuthal  symmetry  in  scattering,  the  first  term  becomes: 

1 


n-0  i 


(10) 


The  integral  is  simply  a moment  equivalent  to  a Legendre  coefficient,  ^ in 

D 

the  flux  expansion: 


♦(»•)  - I (^)*.  p.(ii') 


t^O 


2 ''•'n  n' 


(11) 


10 


Therefore  the  first  part  becomes  simply. 


Air  ' n n 

n“U 


The  second  part  can  be  rewritten  to  contain  the  Integral 


d?'  cose(?-c'). 


Making  the  substitution. 


a-e(C-C’), 

the  Integral  can  be  rewritten  as 


- -p  cos  ada  ■ 4 [slna]  - 0.  (15) 

^ ® 3(C-2ir) 

iKC-2ir) 

Therefore  the  entire  scattering  term  reduces  to  the  simple  form  of  expression 

(12). 

If  It  were  desired  to  express  the  cross  section,  S (Q*-*{2)  in  some 

8 

form  other  than  Legendre  polynomials,  presumably  one  In  which  relatively 
low  order  were  possible,  but  continued  to  express  flux  In  Legendre 
polynomials,  then  the  scattering  term  would  take  the  form: 


E,<5-^)  ^ I r (M')  ia- 

n-0 


1 ” 

<2tr*-lH  E (S'-^)  P„(y’)dn' 


11 


The  moments  defined  by  the  Integral  would  be  different  for  each 
final  direction  considered.  Thus,  for  Instance,  In  a discrete-ordinate 
formulation  using  m quadratures,  a total  of  m times  n moments  would  be 
necessary.  If  this  could  be  shown  to  be  less  than  the  present  n necessary 
for  Legendre  expansion  of  cross  section,  the  situation  might  be  Improved. 

If  a completely  different  orthogonal  series  were  used  to  expand 

both  cross  section  and  flux,  then  an  expression  similar  to  (7)  would  have 

to  be  found  for  the  series  used  or  the  expansion  In  cross  section  carried 

for  all  combinations  of  Initial  and  final  direction.  That  Is,  If  the 

equivalent  of  P (y  ) In  the  new  polynomial  were  to  be  retained,  then 
n o 

Instead  of  using  polynomials  In  terms  of  a single  set  of  directions  y, 
polynomials  would  have  to  be  evaluated  In  terms  of  all  combinations  of 
Initial  and  final  directions. 

2 

If  m were  the  number  of  directions,  m Instead  of  m polynomials 
would  be  used  and  n times  m cross-section  coefficients  evaluated.  As 
before.  If  this  could  be  shown  to  be  less  than  the  n presently  used  In 
Legendre  expansions,  an  Improvement  could  result. 


i • 
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4.2 


Legendre  Coefficient  Data  Required 

The  differential  scattering  cross  section,  a(6)  Is  the  probability 
of  scatter  per  unit  solid  angle  about  6.  Then  the  total  probability  of 
scattering  (of  a unit  light  flux  per  scattering  particle) Is 

IT  IT 

■ I a(9)dfl  ■ 2ir  j a(6)sln6d(t> 

0 0 


1 


« 2ir  j 
-1 


O(0)dx  . 


(22) 


It  Is  desired  that  a (9)  be  expressed  In  terms  of  Legendre  polynomials  and 
a set  of  coefficients: 


a(0)  - [ C P (x) 
n n 


n»0 


(23) 


or,  to  sufficient  approximation, 

N 

a(0)  = l c ?(x)  . 
n-O  ° 


(24) 


These  coefficients  can  be  evaluated  from  a(0)  Information  by 
Integrating  over  appropriate  polynomials: 


0(0)  P„(x)dx  ~ I C 
” n-0  “ 


P (x)P  (x)dx 
n m 


-1 


-1 


(by  the  orthogonality  relationship 
of  the  Legendre  polynomials.) 


(25) 


14 


'-sT', 


c - (^)  a(0)P^(x)dx. 

n * n 


Since  the  n-0  polynomial,  P (x)“l,  from  Eq.  (26), 

o 


Alternatively,  the  set  could  be  defined  by 


o'  = Air  C 
n n 


<=„-|  • 


A set  of  normalized  coefficients  O can  be  defined  by 

n 


a = 2rr  (^r-rr)  C 
n 2n+l'  n 


where  a •*  a . 

O 8 


where 


c(6)-  -L  I p„(.) 

n«0 


o ■ 2Tr  [ o(6)  P (x)dx 
n J n 

-1 


In  which  case 


0(0)  = y o'  P (x) 
' ' 4n  n n 
n»0 


1 

) O(0)P^( 

- 


O^  - (2n+l)a^  - 2Tr(2n+l)  O(0)p^(x)dx. 


This  is  the  case  assumed  In  AMISN  and  other  Oak  Ridge  codes. 


Thus,  where  the  coefficients  are  defined  by  Eq.(28),  they  must  be 
divided  by  (2n-l-l)  for  Input  to  codes  assuming  definition  (31).  Generally, 
Los  Alamos  codes  (DTP,  TRANZIT,  TWOTRAN)  assume  definition  (28). 

4.3  Application  of  Mle  Theory  Legendre  Coefficients  to  ANISN 

In  ANISH-related  terms,  the  fraction  of  scattering  In  a solid 
angle  about  6 Is 


am 


1 ** 

“ I a'  T (x)  . 

4tio_  "_.^n  n'  ' 


s n-0 


This  Is  equivalent  to  f(0)  of  Mle  theory.  From  Eq.(18), 

n«0 


Therefore  the  coefficient  required  by  AMISN,  o',  can  be  determined  by 

n 

equating  the  two  relationships: 


(34) 


(35) 


, N N 

7-  y a P (x)  - 7^—  y o'P  (x) 

4Tr  ^ _ n n ' ^ /I'ltn  “ n n ' ' 


4iTa  ■ n n ' 
8 n**0 


(36) 


The  equality  Is  satisfied  If,  for  all  n. 


Since 


so  that 


o'  ■ a o . 
n ns 


irr 


I 2 * 

a'  - Kmr  a 
n n 


(37) 


(38) 


For  use  In  LASL-type  codes  using  the  DTF  format  (defined  by  Eq.  (30)),  Mle 
theory  coefficients  must  be  divided  by  (20^1).  This  operation  Is  carried 
out  In  the  XSLA  option  of  MACRO  or  MIELE6. 
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The  Mle  Legendre  coefficients,  a , are  a function  of  particle  size;  or 

more  correctly,  of  the  size  parameter,  a - , where  X Is  the  radiation 

wavelength.  Then  for  any  distribution  of  particle  sizes,  n(r),  the 

coefficients  used  In  the  ANISN  transport  calculation  must  use  appropriate 

averages  over  the  distribution.  The  dimensions  of  <J'  are  those  of 

2 ” 

microscopic  cross  section  (distance)  . Ultimately,  ANISN  requires  data 
In  dimensions  of  macroscopic  cross  section,  (distance  ^).  This  is  achieved 
by  multiplying  the  microscopic  cross  section  by  scattering  particle  density. 
If  the  distribution,  n(r),  is  density  per  unit  of  particle  size,  such  that 


N 


■ 

n(r)dr 

I 

0 


(39) 


gives  the  total  particle  density,  then  the  appropriate  macroscopic  coefficients 
that  ANISN  can  use  directly  In  transport  calculations  are 


00 


z - 

n 


o' (a)n(r)dr. 
n 


0 


(40) 


Or,  for  calculation  of  Z where  values  of  a'  are  available  for  specific 

n n 

values  of  a. 


where 


Ar^  Is  the  range  of  sizes  represented  by  r^  and 


(41) 


Expanded  to  angular  cross  section, 

n«0 


(42) 


V.  INFORMATION  FLOW 


In  general,  the  phase  function  coefficients,  a^,  can  be  used 
to  generate  Interaction  probability  coefficients,  A^,  for  different 
combinations  of  size  and  wavelength  .through  the  size  parameter,  a. 
Also,  for  large  values  of  a,  the  computer  time  necessary  to  calculate 
values  of  a^  Is  sufficiently  great  that  unnecessary  recalculations 
should  be  avoided.  Therefore,  the  procedures  for  preparing  Input 
data  for  a discrete  ordinates  code  Is  built  around  a file  of  collected 
values  of  a^.  This  file,  LEG,  can  contain  all  values  previously  cal- 
culated and  can  be  drawn  upon  to  produce  values  of  A^  for  discrete 
ordinates  Input  data.  A diagram  of  the  Information  flow  between  codes 
and  files  Is  shown  In  Fig.  4. 


The  MIELEG  code  Is  based  on  the  MIE-2  code  (Ref.  7)  which  cal- 
culates Mle  Theory  cross  sections  and  phase  functions,  and  the  Legendre 
expansion  expressions  outlined  by  Clark,  Chu,  and  Churchill  (Ref.  8). 
This  code  calculates  from  r,  X,  and  3 a set  of  values  of  a and  adds 
them  to  the  LEG  file.  Any  other  set  of  r,X  values  leading  to  the  same 
size  parameter,  a,  would  produce  the  same  set.  For  a range  of  values 
of  r,  a number  of  sets  can  be  produced,  distributed  across  the  range. 

In  addition  to  saving  these  sets  on  the  LEG  file.  If  the  parameters 
of  a distribution  function  are  Input,  MIELEG  will  produce  a set  of 
macroscopic  parameters,  A^.  These  will  be  calculated  from  the  approxi- 
mation to  Eq.  (4)  : 


A 


n 


NRS  , 

I K(a^,3)  Trr^  a^(oi^,3)  n(r^)  Ar^ 


(43) 


where 


1 

Is 

NRS 

is 

is 

n(r^) 

Is 

the  Index  of  a single  r, 

the  number  of  r values  used, 

the  Interval  of  r associated  with  value  r^, 

the  relative  density  of  particles  of  size  r^. 


These  values  of  A are  written  on  the  MAC  file  for  use  In  making  up 
n 

data  sets  for  discrete  ordinates  transport  calculations.  (Values  are 
written  on  E12.5  format,  so  can  be  used  directly  In  most  codes. 
Including  free-format  ANISN  and  TDA.) 


For  many  A^  sets  of  Interest,  It  Is  unlikely  that  all  contributing 
a^  sets  would  be  calculated  In  a single  MIELEG  run.  It  Is  far  preferable 
to  generate  the  required  macroscopic  A^  set  from  the  collection  of  micro- 
scopic a^  sets  on  the  LEG  file.  This  not  only  Is  more  conservative  of 
computer  resources,  but  should  give  better  approximation  of  Eq.  (43)  to 
Eq.  (4)  by  allowing  a finer  definition  of  r-values.  Macroscopic  sets 
can  then  be  generated  from  the  LEG  file  by  the  MACRO  code,  which  Is 
essentially  the  last  part  of  the  MIELEG  code. 


All  MIELEG  runs  automatically  save  a^  values  on  the  LEG  file. 

New  sets  are  saved  at  the  beginning  of  the  file.  (Old  values  are  copied 
to  a temporary  file  containing  the  new  ones,  then  copied  back  to  the 
original  file.)  Since  MACRO  expects  sets  to  be  ordered  by  ascending 
values  of  the  size  parameter  a,  the  new  LEG  file  must  be  reordered 
after  contributions.  An  ordering  code,  LEGEDIT,  properly  orders  the 
file  and  deletes  old  sets  when  new  ones  with  the  same  value  of  a are 
encountered.  The  LEGEDIT  code  Is  Intended  to  be  used  separately  or  In 
the  same  run  string  with  MIELEG,  thereby  assuring  that  an  unordered  LEG 
file  never  be  encountered. 


In  general,  one  would  need  to  have  a^  values  resulting  from 
different  3-values  on  different  files.  The  LEGEDIT  code  compares  the 
3-value  of  all  sets  with  that  of  the  last  one  on  the  file.  Any  set 
whose  value  does  not  agree  with  that  of  the  last  set  Is  deleted  and 
written  on  an  auxiliary  file,  OTHERLEG.  It  Is  thus  possible,  by 
changing  file  assignments  In  the  LEGEDIT  Job  control  sequence,  to 
unfold  a file  containing  sets  with  several  3-values  Into  several  sets 
each  with  a single  3-value. 


^ In  ordinary  production  of  a set  of  values  for  use  In  a 

transport  code,  one  would  choose  the  LEG  file  containing  the  3-value  i 

of  Interest  and  Input  to  the  MACRO  code  the  wavelength  and  size  dls- 
trlbutlon  required.  MACRO  will  then  produce  a file,  MAC,  containing 
the  set  In  E12.5  format.  This  file  can  be  edited  by  CRT  terminal  to 

MERGE  the  A^  values  Into  the  transport  code  data  file.  | 

Mle  theory  Legendre  coefficient  sets  are  quite  lengthy,  j’ 

requiring  a much  higher  order  than  Is  necessary  for  transport  of  jj 

other  radiations.  Proper  use  of  discrete  ordinates  codes  generally  i 

requires  S-orders  (angular  quadrature)  greater  than  or  equal  to  P-orders 
(Legendre) . (Excessively  large  S-order  may  require  very  long  computing 

times.)  Thus,  for  convenience  In  changing  quadrature  sets  (weights  |' 

[ and  cosines),  a file,  QUADS,  has  been  established.  Into  which  all  Input 

sets  can  be  read  for  merging  Into  transport  code  Input  files.  ! 

A short  normalizing  routine,  MACEDIT,  has  been  Included  for 

r 

data  normalization  In  terms  of  the  first-order  coefficient  A^.  This  j: 

Is  equivalent  to  the  total  scattering  cross  section  and  can  be  used 
: Co  normalize  the  entire  set.  The  normalized  set  (along  with  the 

original  set  from  MAC)  Is  saved  on  the  MAC2  file. 

All  procedures  were  Implemented  on  the  SI®lA-9  computer  at  i 

' Texas  Christian  University,  accessed  by  CRT  terminal.  All  discussion  ! 

Is  directed  toward  usage  on  this  system.  j 


VI.  CROSS-SECTION  ANALYSIS 


A significant  characteristic  of  Legendre  expansions  of  Hie  theory 
cross  sections  is  the  large  order  (l.e.,  number  of  coefficients)  that  is 
necessary  for  adequate  representation.  The  order  needed  is  about  2a,  where 
a is  the  size  parameter,  2TrR/X.  (Some  authors  suggest  2cH-2.)  If  the 
expansion  is  performed  to  produce  angular-dependent  from  order-dependent 
data  (Eq.  1 or  Eq.  3),  a final  term,  N,  greater  than  2oH-2  should  not 
contribute  significantly  to  the  result.  As  examples  of  this  characteristic, 
three  sets  of  microscopic  coefficients,  a^,  produced  by  MIELEG  are  shown  in 
Pig.  5. 

An  indication  of  the  sufficiency  of  the  expansion  can  be  obtained 
by  refitting  the  expanded  function  by  integrating,  for  each  order,  the 
function  multiplied  by  the  Legendre  polynomial  of  that  order  (in  the  manner 
of  Eq.  24) . In  Fig.  6 is  a plot  of  a weighted  average  of  coefficients  from 
o>8.46  to  a ~ 25.22  (normalized  to  1.0)  and  the  coefficients  obtained  from 
a refit  of  the  expanded  function. 

The  refitted  coefficients  become  effectively  constant  when  the 
original  coefficients  become  small.  The  calculation  was  performed  using 
two  different  recursion  relations  for  generation  of  the  Legendre  polynomials, 
and  with  both  single  and  double  precision  in  the  polynomial.  Since  there 
were  no  significant  differences  in  the  resulting  coefficients,  one  must 
assume  that  the  nearly-constant  values  result  from  the  availability  of 
detail  in  the  expanded  function,  rather  than  from  round-off  error. 


Re-expansion  from  the  refitted  coefficients  would  result  in  little 
contribution  from  the  constant  coefficients  since  the  sum  of  several 
polynomials  is  approximately  zero.  If  the  coefficients  were  to  be  obtained 
by  fitting,  as  from  MIE-2  phase  functions,  rather  than  directly  from  Hie 


Theory  (Eq.  17),  this  "tall"  of  constant  values  would  result.  Judicious 
truncation  of  a set  of  coefficients  obtained  by  either  technique  results 
in  adequate  representation  of  the  function. 

Early  truncation  of  the  series  may  give  a greatly  inadequate 
representation  of  the  function.  Figure  7 is  a comparison  of  a full 
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expansion  of  52  average  coefficient  terms  with  an  expansion  using  only 
the  first  29.  It  Is  Interesting  to  note  that  the  number  of  peaks  and 
valleys  In  the  truncated  expansion  Is  equal  to  the  number  of  terms 
used.  (The  valleys  do  not  appear  on  a semi- log  plot.)  As  larger 
numbers  of  terms  are  used,  the  number  of  peaks  and  valleys  Increase 
and  gradually  converge  on  the  complete  expansion. 


The  use  of  a high  Legendre  order  Imposes  on  a transport  code 
larger  computer  storage  and  running  times  than  one  with  low  order.  It 
Is  likely  that  the  detail  of  the  large  forward  scattering  peak  (the 
phase  function  near  zero  degrees)  more  strongly  affects  transport 
calculation  results  than  the  detail  of  the  smaller  peaks.  A smoothing 
option  was  added  to  the  MIELEG  and  MACRO  codes  to  determine  whether 
some  of  the  non-forward  detail  might  be  dispensed  with.  The  refitting 
calculation  discussed  earlier  performed  the  Integration  over  function  and 
polynomial  by  summing  one-degree  Increments.  The  smoothing  option  performed 
the  same  Integration  using  one-degree  Increments  In  the  forward  peak  but 
larger  Increments  at  other  angles.  A macroscopic  set  of  coefficients  was 
determined  for  a range  of  a of  0.85  to  8.5.  This  was  then  e3q>anded  and 
refitted  with  non-forward  Increments  of  10  and  20  degrees.  The  resulting 
coefficients  were  then  re-expanded  to  test  the  representation  resulting 
from  the  "smoothed"  fits.  The  expanded  functions  (macroscopic  cross  sections) 
are  shown  In  Fig.  8. 

Using  10-degree  smoothing,  the  Increments  are  narrower  than  the  peaks 
appearing  In  the  function  and  re-expandlng  the  coefficients  gives  results 
comparable  to  the  original  function.  Using  20-degree  smoothing,  however, 
the  Increments  are  comparable  to  the  peak  width  and  expansion  of  the 
resulting  coefficients  shows  a result  similar  to  that  observed  with  a 
truncated  expansion.  The  function  resulting  from  smoothed  fitting  seems 
to  require  higher  order,  rather  than  less.  Consideration  of  the  method  by 
which  the  coefficients  are  obtained  — Integration  of  the  function  over  all 
angles,  weighted  by  Legendre  polynomials — should  lead  one  to  recognize 
that  the  predominant  effect  must  result  from  the  large  forward- scat ter  peak. 
Thus  the  characteristic  periodicity  of  the  function  Is  determined  by  this 
predominant  peak,  and  Its  damping  for  higher  angles  Is  determined  by  the 
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detail  employed.  In  the  same  way  that  the  order  la  related  to  the  number  of 
peaks,  the  order  must  be  related  to  the  narrowness  of  the  forward  peak. 
Adequate  representation  of  the  forward  peak,  therefore,  requires  a high 
order  and  no  "short-cut"  Is  available  In  Legendre  representations  of  highly 
forward  scattering. 
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VII.  TEST  EXPERIENCE  - C(»<PUTER  REQUIRQIENTS 


7.1  Cross-Section  Generation 

A series  of  cross-section  sets  was  produced  and  tested  on  the  Time 
Dependent  ANISN  (TDA)  code.  A case  was  modeled  with  a unit  point  source  at 
the  center  of  a 0.5-km  radius  cloud.  Time  Increments  summing  to  more  than 
10  ysec  were  specified.  (The  unscattered  wave  front  would  reach  the  edge 
of  the  cloud  In  1.67  ysec) . Early  attempts  to  generate  complete  cross- 
section  sets,  using  particle  sizes  as  would  be  expected  In  a cloud,  required 
long  MIELEG  running  times  (considerably  more  than  the  few  minutes  estimated) . 
In  addition,  the  large  orders  Involved  appeared  to  require  large  computer 
storage  and  CPU  time  with  TDA.  At  this  time  It  was  determined  to  produce 
the  LEG  file  and  LEGEDIT  and  MACRO  codes  to  conserve  cross-section 
preparation  time  and  to  perform  test  calculations  on  ranges  of  particle 
sizes. 


The  particle  density  function  used  In  the  tests  was  given  by  the 
expression 

2 -0.6R 

n(R)  - 27R  e (44) 

where  R Is  given  In  microns  (10  ^ cm)  and  n(R)  has  units  of  particles 
cm  . The  density  function  Is  shown  In  Fig.  9. 

Microscopic  coefficients  were  produced  by  MIELEG  for  47  size 
parameters  up  to  a value  of  84,  and  saved  on  the  LEG  file.  All  assumed  a 
complex  refractive  Index  of  (1.33,  0.0).  For  a light  wavelength  of  0.75 
microns,  the  LEG  file  can  provide  data  for  particle  sizes  to  10  microns. 

The  CPU  time  as  a function  of  size  parameter,  a.  Is  shown  on  Fig.  10.  All 
times  Include  a following  LEGEDIT  run  to  order  the  file. 


7.2  Angular  Quadrature  Order 

The  cross-section  order  (P-order)  necessary  to  reliably  represent 
Mle  phase  functions  was  determined  earlier.  A related  parameter  Is  the 
angular  quadrature  order  (S-order)  necessary  to  adequately  represent  the 
anisotropy  of  scattering  In  the  balance  equations  of  the  discrete  ordinates 


calculations.  The  "quadratures"  are  the  vector  direction  cosines,  x, 
along  which  angular  fluxes  are  calculated.  Associated  with  the  cosines 
are  weights,  w,  representing  the  fractional  solid  angles  around  each 
vector  direction.  Taken  together,  these  represent  discrete  directions 
(thus  "discrete  ordinates")  that  the  radiation  can  take.  Highly  peaked 
cross  sections  generally  require  the  use  of  a large  number  of  directions 
(high  S-order)  to  avoid  loss  of  Information  In  the  calculation.  A general 
rule  of  thumb  Is  that  the  S-order  should  be  no  less  than  the  P-order. 

Several  TDA  calculations  were  performed  with  particle  sizes  of 
0-2  microns  and  a wavelength  of  0.75 microns  (P-order  of  35  for  2ori-2 
representation).  S-orders  of  16,  24,  32,  and  48  were  used.  For  the  sake 
of  comparison,  the  results  with  S-48  were  assumed  to  be  exact.  Fluxes  at 
the  outer  edge  of  the  cloud  (500  m)  were  compared  for  all  time  steps.  Average 
percent  errors  were  calculated  for  1-4  and  4-10  transit  times  (1  transit 
time:  1.67  sec,  time  for  unscattered  wave  front  to  reach  the  edge  of  the 
cloud).  The  results  are  shown  In  Fig.  11.  All  errors  were  negative. 
Indicating  that  the  fluxes  at  the  outer  edge  were  reduced  by  Inadequate 
angular  representation  of  the  radiation  transport.  The  error  was  greater 
for  the  longer  times  considered,  and  appeared  to  Increase  sharply  below 
S-24.  It  appears  that  the  "S-order  greater  than  P-order"  rule  might  be 
broken  for  moderate  error  ('V'5%)  but  could  lead  to  much  larger  error  If 
S-orders  are  allowed  to  fall  much  below  the  P-order. 

All  quadrature  sets  reported  here  were  symmetric,  having  a negative 
cosine  to  match  each  positive  cosine.  It  Is  possible,  exercising  great 
care,  to  i;ise  "asymmetric  quadrature  sets",  which  have  a greater  density  of 
directions  considered  near  x>1.0  to  better  represent  the  forward  scattering. 

A single  TDA  run  was  made  using  such  a quadrature  set  (of  unknown  origin) . 
These  results  are  not  reported  except  to  say  that  the  fluxes  were  20-50% 
greater  than  the  S-48  fluxes,  so  should  not  be  accepted  as  reliable  without 
further  study. 

7.3  Computer  Storage  Requirements 

Several  other  TDA  Jobs  were  run,  with  maximum  particle  sizes  up  to 
6 microns,  with  P-orders  up  to  95  and  S-orders  from  16  to  48.  Core  storage 
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requirements  are  directly  related  to  the  S-  and  P-orders.  The  MAIN 
routine  of  TDA  contains  a COMMON  statement  specifying  the  size  of  the 
flexibly-dimensioned  storage  array.  Load  modules  containing  array  sizes 
of  20,000  (file  name;  TDA)  and  30,000  (TDA30K)  were  used.  It  is  likely 
that  the  largest  array  size  that  can  be  used  with  the  TDA  code  and  the 
SIQIA  9 computer  Is  about  46,000  with  a total  core  usage  of  64,000. 

Any  change  in  the  COMMON  array  must  be  matched  by  the  value  of  the  size 
test  variable  LIMl.  TDA  prints  the  size  of  COMMON  array  necessary  for 
any  problem  at  read-ln  time,  so  the  requirements  can  be  determined  with- 
out actually  running  a job  by  submitting  a job  with  a very  short  run  time 
specification  (<1  minute).  For  the  TDA  test  case,  C0MhK)N  array  require- 
ments for  all  P-S  combinations  tested  are  shown  on  Fig.  12.  (The  P-values 
shown  on  Fig.  12  are  the  values  for  data  stored.  The  P-value%  of  the 
calculation  were  often  less  as  truncation  to  2a+2  was  employed.) 

Using  the  estimated  errors  In  Fig.  11,  It  Is  possible  to  supply  lines 
of  estimated  error  on  Fig.  12,  using  as  the  estimate  for  S-order  equal  to 
P-order,  and  20%  for  S-order  equal  P-order  /2.  Extrapolating  these  results 
to  an  8-mlcron  truncated  P-order  of  >130,  with  an  S-order  of  64,  a storage 
requirement  of  '^2K  would  result  for  20%  error,  near  the  maximum  of  46K  for 
the  Slgma-9  computer. 

The  other  computer  requirement  of  consideration  is  the  running  time 
(CPU  time).  From  the  same  test  runs  the  times  are  shown  on  Fig.  13,  again 
cross-plotted  with  estimated  5%  and  20%  error.  In  this  case,  the  values 
were  plotted  against  the  P-order  actually  used  In  the  calculation. 
Extrapolation  to  a meixlmum  particle  size  of  8 microns  Is  more  difficult  for 
CPU  time  than  It  was  for  COMMON  storage,  but  It  appears  that  a run  time 
of  more  than  100  minutes  for  20%  error  is  likely. 


li 


Successful  use  of  asymmetric  quadrature  sets  could  reduce  the 
effective  S-order  by  a factor  of  3 or  more  while  retaining  the  degree  of 
error.  This  could  reduce  the  storage  requirements  to  '^'1/2  and  time  to 
'^'1/3  for  the  test  cases  reported  here.  It  Is  likely,  then,  that  use  of 
asymmetric  quadrature  sets  could  extend  the  range  of  solution  to  12  microns 
or  more,  effectively  covering  the  particle  size  range. 
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LIKELY  MAXIMUM 


SIGMA-9  COMPUTER 


STORAGE  REQUIREMENTS  FOR 
TDA  TEST  PROBLEM 


Constant  S-Order 


P-ORDER  STORED 


Fig.  12.  COMMON  Array  Storage  Requirements  for  TDA  Test  Problem 


P-ORDER  USED 


Fig.  13.  Running  Time  Requirements  for  TDA  Test  Problem 
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VIII.  TEST  EXPERIENCE  - CALCULATIONAL  RESULTS 


Test  calculations  were  performed  with  the  TDA  code  for  0-1 
micron,  0-2  micron,  0-4  micron,  and  0^-6  micron  particle  ranges  with 
a wavelength  of  0.75  microns,  complex  refractive  Index  of  (1.33,0.0) 
and  a particle  size  distribution  defined  by  the  expression: 

n(R)  - 27R^  (47) 

where 

R Is  the  particle  radius,  microns 

-4 

n(R)  Is  the  density  distribution  In  particles  cm  . 

Macroscopic  Legendre  coefficients  were  generated  for  each 
range  from  Individual  particle  size  parameter  data  stored  on  the  LEG 
file.  The  MIELEG  code  was  used  to  generate  the  Individual  particle 
data  and  the  MACRO  code  to  process  It  Into  macroscopic  data.  Each 
Increase  In  particle  size  range  requires  an  Increase  In  the  Legendre 
order  for  adequate  representation  of  the  macroscopic  function  (phase 
function  Integrated  over  particle  size  distribution).  Also,  as  more 
particles  are  considered,  the  total  cross  section*  (the  first  Legendre 
coefficient)  Increases.  The  necessary  order  and  total  cross  section 
are  shown  In  Table  I. 


TABLE  I.  PARAMETERS  OF  MACROSCOPIC  DATA  USED  IN  TEST 
CALCULATIONS 


Particle  Size 
Range,  microns 

Maximum  Order 

Req'd  (2ori-2) 

Total  Cross 
Section,  cm~l 

0-1 

18 

1.72X10“® 

0-2 

35 

5.59X10“® 

0-4 

70 

5.84X10“^ 

0-6 

103 

1.72X10“^ 

0-8 

136 

3.00X10“^ 

Actually,  this  Is  the  total  "scattering"  cross  section.  If  a non-zero 
extinction  coefficient  were  used,  the  particle  density  times  extinction 
coefficient  would  be  the  "absorption"  cross  section  which,  when  added  to 
the  total  scattering  cross  section,  would  equal  the  total  cross  section. 
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While  orders  of  18,  35,  emd  70  were  used  for  the  1-,  2-,  and  4-inlcron 
calculations,  microscopic  particle  size  data  up  to  only  5.5  microns 
(but  Integrated  to  6 microns)  was  used  to  estimate  the  cross  sections 
to  6 microns,  thereby  allowing  an  order  of  only  95  to  be  used.  The 
8-mlcron  calculation  was  not  performed,  due  to  computer  limitations, 
but  was  Included  In  the  table  to  note  that  that  rate  of  Increase  of 
cross  section  decreases  when  the  range  of  sizes  extends  well  beyond  the 
density  distribution  peak. 

The  number  of  particle  sizes  considered  Influences  the  values 
of  the  macroscopic  coefficients  to  some  extent.  As  a test  of  this, 
the  coefficients  of  the  0-4  micron  range  were  generated  with  the  MACRO 
code  from  LEG  file  data,  using  25  particle  sizes,  and  directly  with  the 
MIELEG  code  using  four  sizes  (PARTS  - 3.0).  The  more  approximate 
MIELEG  results  Included  a total  cross  section  of  6.75X10  ^cm  compared 
to  the  MACRO  result  of  5.84X10  The  coefficients  obtained  In  these 
cases  are  shown  In  Fig.  14.  Whereas  the  MIELEG-calculated  coefficients 
are  greater  than  the  MACRO-calculated  coefficients.  It  Is  expected  that 
some  other  choice  of  the  four  sizes  used  by  MIELEG  might  have  resulted 
In  smaller  coefficients. 

It  was  evident  In  Section  VII  that  the  number  of  significant 
coefficients  used  In  an  expansion  was  approximately  equal  to  the  number 
of  peaks  + valleys  In  the  expanded  function  (cross  section  versus  angle) . 
Since  the  peaks  and  valleys  for  any  single  particle  size  are  located 
at  different  angles  from  those  for  any  other  particle  size,  the  summa- 
tion for  several  particle  sizes  will  result  In  a larger  number  of  peaks 
and  valleys  than  there  would  be  for  the  largest  particle  size  con- 
sidered. Since  the  order  of  the  summation  Is  equal  to  the  order  of 
the  largest  particle  size,  the  function  detail  will  not  be  completely 
represented  by  the  summation.  Evidence  of  this  Insufficiency  may  be 
an  Irregular  character  of  the  coefficient  distribution,  as  seen  In 
Fig.  14  for  the  MACRO-generated  coefficients  In  the  40-60  order  range. 

The  expanded  function  resulting  from  the  MACRO-generated  coefficients 
showed  some  of  the  oscillations  characteristic  of  underdetermination 


seen  for  Insufficient  order  (Fig.  7)  or  loss  of  detail  by  smoothing 
(Fig.  8).  The  oscillations  were  not  so  apparent  In  the  expanded 
function  resulting  from  MIELEG-generated  coefficients,  where  fewer 
particle  sizes  were  considered. 

It  Is  not  necessarily  detrimental  to  a reasonable  transport 
calculation  to  use  coefficients  which  would  lead  to  oscillatory 
functions.  The  oscillations  occur  at  angles  for  which  the  scattering 
Is  relatively  Improbable.  (For  0-4  micron  particles,  the  oscillations 
occur  at  'v<10  of  the  forward  peak  value.)  After  a few  scatters, 
nearly  forward  scattering  would  mask  any  error  resulting  from  the 
oscillatory  behavior. 

TDA  calculations  were  performed  on  the  test  case  using  both 
sets  of  coefficients.  The  fluxes  obtained  at  the  edge  of  the  0.5-km 
radius  cloud  from  a point  source  at  the  center  emitting  a unit  pulse 
at  time  zero  are  shown  In  Fig.  15.  Also  shown  Is  the  case  for  which 
the  MIELEG  coefficients  are  normalized  so  that  the  total  cross  section 
(first  coefficient)  Is  the  same  as  the  total  cross  section  from  MACRO. 
Except  for  shortest  times  ('^>1  ysec  after  arrival  of  the  unscattered 
wave  front  at  the  edge  of  the  cloud) , the  normalized  result  falls 
between  the  MIELEG-  and  MACRO-based  results.  This  Indicates  that,  for 
this  case,  the  difference  In  values  of  total  cross-section  and  the 
difference  between  the  angular  details  are  about  equally  Important  In 
determining  the  flux  at  the  edge  of  the  cloud.  It  Is  not  known  whether 
the  angular  detail  Is  better  for  the  MIELEG  data,  based  on  few  particle 
sizes,  or  for  the  MACRO  data,  which  appears  underdetermined  for  low 
probability  scattering.  The  value  of  the  peak  flux  Is  not  reliable 
as  Insufficient  detail  was  used  In  spatial  Increments  at  the  outer 
edge  and  time  Increments  when  the  wave  front  reaches  the  outer  edge 
to  provide  good  definition. 

The  flux  at  the  outer  edge  of  the  cloud  Is  plotted  versus  the 
time  after  the  unit  pulse  In  Fig.  16  for  all  particle  ranges  considered. 
All  results  shown  here  were  calculated  with  2ori'2  truncation  of  P-order 


fluxes  at  edge  of  0.5  KM  CLOUD 
0-4  MICRON  PARTICLES 


MIELEG  Data  (4  sizes) 


MACRO  Data 


TIME  AFTER  PULSE,  usee 


Fig.  15.  Effect  of  Coefficient  Data  Source  on  TDA-Calculated  Fluxes 


and  S48  quadratures.  The  error  in  the  results  from  using  Insufficient 
S-order  Is  likely  to  be  about  20Z  for  the  0-6  micron  case  at  longest 
times  and  less  for  all  other  times  and  size  ranges.  It  Is  clear  that 
a greater  range  needs  to  be  considered  to  observe  the  effect  of  the 
entire  range  of  sizes  In  the  distribution. 

With  the  computer  limitations  of  this  st\idy,  extension  to 
sufficiently  greater  sizes  would  require  a significant  decrease  In 
P-order  or  S-order.  The  possibility  of  decreasing  the  S-order  by  the 
use  of  asymmetric  quadrature  sets  has  already  been  mentioned  (Section 
VII).  The  separation  of  total  cross  section  and  angular  detail,  as 
practiced  In  normalizing  MIELEG-calculated  coefficients  to  the  total 
cross  section  from  MACRO,  provides  a possibility  for  modest  reduction 
In  P-order. 

For  transport  cases  of  sufficient  depth  so  that  the  forward 
direction  of  the  Initial  pulse  becomes  greatly  spread,  a precise 
definition  of  forward  scattering  might  not  be  necessary  for  adequate 
results.  Less  forward  scattering,  as  represented  by  lower  order  repre- 
sentation, would  predict  a longer  path  to  locations  of  interest,  so 
lower  order  results  might  be  corrected  by  multiplying  the  time  scale 
by  a factor  equal  to  the  ratio  of  expected  path  lengths.  As  a brief 
teat  of  this  possibility,  two  low-order  angular  distributions  of 
cross  section  were  normalized  to  the  total  cross  sections  of  two  other 
cases  of  higher  order.  The  chosen  cases  are  shown  In  Table  II. 

TABLE  II.  TEST  CASES  FOR  LOW-ORDER  REPRESENTATION  OF 
HIGH-OROER  SCATTER 


Low-Order  Cases 

Hlgh-Order  Cases 

Rayleigh  scatter 

0-4y  (MACRO) 

0-2y  (MACRO) 

0-6u  (MACRO) 
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Rayleigh  scattering  represents  a low-order  limit  of  light 
scattering  and  Is  given  simply  by  the  expression: 

0(6) /Os  " xlif  (1  + cos^e)  (48) 

where  a(6)  Is  the  angular  microscopic  cross  section  and  a Is  the  total 

s 

microscopic  cross  section.  Their  ratio  Is  the  phase  function.  A 
Legendre  representation  of  this  Is  obtained  by  Integrating  It  over  all 
solid  angles,  weighted  by  the  Legendre  polynomial  of  order  n,  to  obtain 
the  coefficient  of  order  n.  Thus, 


• Ztt 


1 

■ 

P (cos0)d(cos0) 

-1 


(49) 


The  first  three  Legendre  polynomials  are  (using  computer 
indices:  1,  2,  3 for  the  first  three  pol3momlals) : 


Pj^(cos0)  » 1 


P2(co80)  “ cosS 


P2(co80)  - I cos^O  - i 


(50) 


Then  a^^  ■ 1,  a^  ■ 0,  a^  - j.  All  further  terms  are  zero.  That  this 
Is  an  exact  representation  is  confirmed  by  substituting  the  coefficients 
Into  the  expansion: 


(51) 


_1_ 

4it 


1.1+0  • CO80  + 


ll  - 2) 
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the  original  expression.  A Rayleigh  scattering  angular  distribution 
can,  therefore,  be  used  to  represent  a higher  order  case  by  a P3  data 
(P2  In  traditional  Indexing) , where  the  first  coefficient  Is  the  total 
macroscopic  cross  section,  the  second  Is  zero,  and  the  third  Is  one- 
half  of  the  total  cross  section.  This  was  done  for  both  0-4  micron 
and  0-6  micron  cases.  In  addition,  the  0-2  micron  P35  angular  distri- 
bution was  normalized  to  both  the  0-4  micron  and  0-6  micron  cases, 
using  the  MACEDIT  code.  The  results  of  both  representations  of  the 
0-4  micron  case,  along  with  the  original  0-4  micron  results,  are  shown 
In  Fig.  17.  Similarly,  results  for  the  0-6  micron  case  are  shown  In 
Fig.  18. 


In  Fig.  17,  both  normalized  cases  have  Initial  peaks,  probably 
from  local  scatter  of  the  passing  wave  front,  followed  by  multl-scatter 
contributions.  The  Rayleigh  scatter  case,  however.  Is  depressed  for 
times  Immediately  after  the  passage  of  the  front,  lacking  a large  com- 
ponent of  strongly  forward-scattered  radiation.  Thereafter,  It  has  a 
very  strong  delayed  component  from  other  scatters.  The  case  normalized 
from  the  0-2  micron  distribution  Is  Intermediate  to  the  two  others, 
having  more  forward-scattering  contribution  than  the  Raylalgh  and  less 
than  the  true  0-4  micron  case.  Within  about  4 ysec  of  the  pulse,  the 
flux  then  appears  of  the  form  of  the  true  case,  but  delayed  by  longer 
scattering  paths. 


In  Fig.  18,  the  Rayleigh  scatter  case  shows  no  Initial  peak  at 
all.  Indicating  that  the  cross  section  Is  large  enough  so  that  the  wave 
front  Is  very  weak  and  local  scattering  from  It  near  the  edge  of  the 
cloud  Is  small.  All  Rayleigh  scatter  contributions  are  delayed  by 
relatively  long  multiple-scattering  paths.  The  case  normalized  from 
0-2  microns  Is  similar  In  shape  to  the  true  0-6  micron  case,  but 
delayed.  This  shows  some  promise  that  the  technique  of  using  a nor- 
malized lower-order  cross  section  set  and  an  adjusted  time  scale  might 
give  reasonable  approximation  to  the  true  higher-order  case  for  large 
numbers  of  scatters. 


cm  sec 


FLUXES  AT  EDGE  OF  0.5  KM  CLOUD 
0-6  MICRfflJ  PARTICLES 


Normalized  from 
0.2  Micron  Scatter 


Calculated 
from  0-6u 
data  (MACRO) 
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APPENDIX  A.  MIELEG  UTILIZATION  INSTRUCTIONS 


The  MIELEG  code  calculates  Legendre  coefficients  of  Interaction  j 

probability  per  unit  particle  density  (microscopic  cross-section 

coefficients,  a^)  or  Interaction  probability  per  unit  path  (macro-  j 

scoplc  cross-section  coefficients,  A ).  Input  data  are  particle  sire  ] 

1 

range  and  distribution,  radiation  wavelength,  and  complex  refractive  ’ 

Index.  It  Is  based  on  the  MIE-2  code  and  direct  evaluation  of  Legendre  ! 

coefficients  (Refs.  7 and  8).  It  produces  a file,  LEG,  of  microscopic 
coefficients  and  a file,  MAC,  of  macroscopic  coefficients. 

A.l  MIELEG  Input  Description 

1 

RMIN,  RMAX,  PARTS,  PM  (Complex),  WAVEL,  NEQ  (6F10.0,  15)  j 

RMIN  - Minimum  particle  radius  considered,  microns  I 

RMAX  - Maximum  particle  radius  considered,  microns  i 

PARTS  - Number  of  Intervals  In  range  of  considered  radii. 

If  0.0  used,  only  RMAX  considered,  no  size  dis- 
tribution necessary,  no  MAC  file  produced. 

PM  - Complex  refractive  index,  8,  composed  of  two  parts:  ] 

a)  M,  real  Index  of  refraction  relative  to  the  j 

surrounding  material,  and 

b)  K,  extinction  coefficient  of  the  particle  material.  j 

WAVEL  - Wavelength  of  the  radiation,  microns  | 

NEQ  - Number  of  equation  describing  particle  size  dls-  I 

trlbutlon,  0 or  1 j 

0:  N(R)  - P1*R**P2*EXP(-P3*R**P4) 

1:  N(R)  - PI  for  R < P2  j 

N(R)  - P3*R**P4  for  R > P2  | 

PI,  P2,  P3,  P4  (4F10.0)  (Ignored  if  PARTS  - 0.0)  j 

Parameters  described  by  choice  of  NEQ,  on  previous  Input  record.  I 


Oitput  Options  4(1X,A4)  (Ignored  if  PARTS  ■ 0.0) 
Possible  options: 

XSOR  or  XSLA,  TRUN,  AGIN,  SMOO 


Output  options  are  entered  In  Hollerith  (alphabetic) form.  They 
can  be  In  any  order,  but  must  be  spaced  as  Indicated  by  the  IX. 
They  can  be  conveniently  Input  with  a leading  space,  then 
separated  by  commas . 

XSOR  or  XSLA  - MAC  option.  Most  codes  requiring  Legendre 
coefficient  cross-section  data  are  of  the  Oak  Ridge  (FIDO) 
type  or  Los  Alamos  (DTP)  type.  The  two  generally  differ  by  a 
factor  of  (2N+1)  where  N Is  the  order  (Initial  order:  zero) 
of  the  Nflst  coefficient.  For  a discussion  of  the  evaluation 
of  coefficients  from  Mle  theory  and  the  type  choice,  see  Section 
IV:  Conversion  of  Mle  Parameters  to  Legendre  Coefficient  Data. 


Option:  XSOR  XSLA 

Type:  Oak  Ridge  Los  Alamos 

Typical 

Codes:  ANISN,  TDA  DTF 

DOT  TWOTRAN 

MORSE  TRANZIT 


TRUN  - Truncation  option. 

The  MIELEG  code  produces  coefficients  for  a higher  order  (more 
values)  than  Is  normally  necessary.  The  MIELEG  code  uses  all 
of  these  coefficients  to  calculate  the  macroscopic  cross- 
section  of  equation  (3)  (labeled  F on  the  printed  output)  for 
0 through  180  degrees.  It  then  truncates  the  series  of  coeffi- 
cients imd  recalculates  the  cross-section  values.  Truncations 
are  chosen  according  to  several  criteria: 

1)  where  the  coefficients  Is  first  smaller  than  the 
Initial  coefficient 

2)  where  N Is  first  greater  than  20+2 

3)  at  local  minima  In  the  coefficient  sequence. 


The  cross  sections  resulting  from  each  of  these  truncations 

can  be  Inspected  to  find  an  adequate  order.  Truncations  that  1 

are  too  severe  often  result  In  a few  negative  cross-section  j 

values.  Truncation  at  2ori-2  Is  generally  adequate. 

AGIN  - Refitting  option.  i 

Test  of  the  consistency  of  the  Legendre  polynomials  and  I* 

I 

expansions.  This  option  Integrated  the  expanded  cross  sections 

over  all  angles,  weighted  by  Legendre  polynomials  to  recalcu-  L 

late  the  coefficients.  In  general,  the  values  agree  with  the 

i 

originals  to  about  3-place  accuracy  until  N^2a,  after  which  the 
original  coefficients  have  contributed  little  to  the  cross 

sections,  so  cannot  be  extracted  from  them.  i 

SHOP  - Smoothing  option. 

This  option  was  added  to  refit  cross-section  data  using  smoothed 
values  for  data  at  angles  greater  than  those  In  the  large 
forward-scattering  peak.  Tests  with  this  option  did  not  result 
In  significantly  different  coefficient  sets.  Indicating  that 
the  Legendre  order  Is  determined  by  the  forward-scatter  peak . 

THRESHF,  WIDTH  (2E10.0) 

Smoothing  parameters.  If  the  SHOP  option  Is  used,  those  para- 
meters are  used.  Cross  sections  are  averaged  over  WIDTH  number 
of  degrees  for  all  angles  greater  than  the  angle  for  which  the 
cross-section  value  falls  to  THRESH*F(1)  where  F(l)  Is  the 
cross-section  value  at  zero  degrees  (the  maximum  of  the  forward 
peak) . 

A. 2 MIELEG  Putput  Description  - Printout 

The  MIELEG  code  prints  out  all  Input  data,  then  radii  and  rela- 
tive densities  If  more  than  one  calculated.  For  each  radius,  the  size 
parameter  (ALPHA),  cross-section  ratio  K(a,3)  (SCKRP) , and  Legendre 
order  (LESS)  are  printed,  followed  by  the  actual  Legendre  coefficients. 


If  several  sizes  are  considered,  the  macroscopic  coefficients  are 
printed  and  these  are  expanded  Into  the  angular  dependent  cross  sections, 
which  are  also  printed. 

If  the  TRUN  option  Is  used,  the  expansion  Is  printed  for  each 
truncation.  If  the  SMOO  and  AGIN  options  are  used,  the  resulting 
coefficients  are  printed  along  with  the  full  expansion.  The  TRUN  option, 
along  with  SMOO  or  AGIN  also  causes  the  truncated  expansions  for  these 
cases  to  be  printed. 

A. 3 MIELE6  Output  Description  - Files 

The  MIELEG  code  produces  two  files,  LEG,  containing  coefficients 
for  single  particle  sizes,  and  MAC,  with  coefficients  Integrated  with 
particle  size  distribution.  Other  files  may  be  produced  by  a MIELEG 
run.  In  case  a run  Is  Interrupted,  the  LEG  file  existing  before  the 
run  would  be  found  on  BACKLEG.  The  new,  partially  completed  file, 
might  be  found  on  NEWLEG.  If  the  run  string  were  completed,  and  Included 
a run  of  the  LEGEDIT  code,  any  coefficients  created  which  had  a different 
complex  refractive  Index  from  those  already  on  the  LEG  file  would  be 
written  on  the  OTHERLEG  file. 

A. 4 MIELEG  FORTRAN  Listing 


ooo  nooooonnnr>  o r>  oo  o ->  nnno 


ROUTINE  MI6LEG  - NEW  {1/7/77)  RECURSION  SCHEME 

COnE  TO  FIT  MIE  phase  FUNCTIONS  WITH  LEGENDRE  POLYNOMIALS 
complex  A(500) .B(500) ,PSIR(500) .PSIRC(500) .BCA(2) .ACAPB, 

1 8ETAI .PM.PSRC.ZETA 

DIMENSION  PSl (3) ,CHI (3) .COEF(500) » 

j R1 (pOO) ,ENR(200) 

COMMON  alpha, Pl  I 

COMmON/POuT/RmIN,RMAX,Pm*WAVEl*NEO,P1 ,P2»P3,P4.NRS, 

1 LESS.ENTOT ,COEFXS(500) 

DIMENSION  SUMW(500) .SUMV(500) 

EOUIVALENcEIA(1).PS1R(1)  ) , «B(1)  .PSIRCd  ) ) 

G(J,K)  = (REAL(A(K)  )*REAL(A(  J)  ) + A IMAG  ( A ( K ) ) * A IMAG  ( A ( J ) ) -t- 
IREAI  (P(K))*REAHB(J))  + A I M AG  ( B ( K ) ) * A I MAG  ( B ( J ) ) ) 

P »<2.*Jd.)/(J»(J  + l. ))*(?.  *K  + l.)/fKa(K  + l.)) 

F(J.K)  =(REAL(A(K> )*REAL(B( J) ) + AIMAG( A(K) )«AIMAG{B< J)  ) + 
1REAL(B(K)  )*REaI  (A(J))  + A IM AG ( B ( K ) ) * A I M AG ( A ( J ) ) ) 

P afP.flJ>l.)/(ja(J+l.))ff(2.*K+l.)/(K»(K+l.)) 

C(K,J)  r (2*N  + 1 )«  ( J*(  J-d  ).fKe(K+l  )-N*(N  + l ) )»»2/4  . 

REWIND  11 
REWIND  1? 

PI  I = 4.0  • ATANC 1 . 0 ) 

PI?  = 2.0  • PM 

DO  2 1=1,500 
CnEFXSM)  = 0.0 
2 'ONTINUE 

read  (5.600)  RM1N,  RMAX,  PARTS.  PM,  WAVEL.  NEO 
WRITE(6.6lO)  RMIN,  RMAX,  PARTS.  PM,  WAVEL,  NEQ 
NRS  - PARTS  * l.O 
IF  (PARTS. OT. 0.0)  GO  TO  25 
RI (1)=RMAX 
GO  TO  72 

USE  PARTS. EO. 0.0  TO  GET  COEFS  FOR  SINGLE  SIZE  ONLY. 
truncating,  refitting  TESTS  SUPPRESSED.  RMIN  IGNORED. 

NO  SIZE  DISTRIBUTION  DATA  REQUIRED. 


25  DELR  = ?.0*(RMAX-RMIN)/(PARTS»(PARTS  + 1 .0)  ) 

*•*  CALCULATE  PARTICLE  RADII  *•* 

Rid)  = RMIN 

DO  28  I=?,NRS 
II  = I - 1 

-’8  RI  ( I)  = RK  ID  ♦ II  aPELR 

*•*  calculate  particle  size  distribution  ••• 

IF (NFO)20.32,40 
<2  READ(5-620)  P1.P?.P3,P4 
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onn  r>  n noon  ooo 


WRIT(:(6.630)  Pl.P2.P3.P4 
no  37  1=1. NRS 

< / ENR( I ) =P1 •RT ( 1 ) «eP?*EXP(-P3*RI { I )#«P4 ) 

GO  TO  70 

40  RFAD(5.620)  P1.P?.P3,P4 
WRl TE  r* .640  )Pt. P2 .P3.p4 

C 

DO  18  1=1, NRS 

IFCRT (!  ).GT .P2)  GO  TO  57 
FnR(I  ) = PI 
GO  TO  18 

57  EnR ( \ ) =P3*R1 ( I )»*P4 
'8  continue 
C 

70  WPITE(6.650)  NRS 

WPI TE( 6,661 ) (RI n ) .ENR( I ) . 1=3 ,NRS) 

C 

/2  ENT0T=0.0 
TEXT  =0.0 
SOOF  =0.0 
C 

no  500  KCsl.NRS 

«**  calculate  mje  variables  *»* 

EXKRO  = 0.0 
SCKRO  = 0.0 
ALPH  = WAVEL/PI2 

71  X = Rl  (KO/ALPH 

TEST  TO  AVOID  SINGUIARITv  AT  MULTIPLES  OF  PI 
IF( A0S<S1N(X) 3 .GT .0.001  ) GO  TO  73 
RI  (KC)=RI {KC)*1.003 
WRITF(6. 74  )KC 

74  FORMAT(/.»R  NUMRER',13.’  FUDGED  TO  AVOID  SINGULARITY  • 
1 'AT  mult.  OF  PI  • ) 

Go  TO  /I 

change  in  X TO  AVOID  ODD  RESULT  at  SINGULARITY. 

73  XYZ  = X 
ALPHA=X 

WOP  : 3 . OE-08  *ALPH»*? 
roFPH=2./(X»X) 

IFlx.l  T.IOO. ) GO  TO  10 
I =X«CABS(PM) 

K1  =1.1*L 
GO  To  11 

10  L rl5. 0*SORT <1 .♦O. 01 *X#X) 

K1  = t .5  • L 

11  K?  = Kl-1 
IF(L.LT.501)  go  TO  1? 
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c 


c 


WRITF(6,670  ) L tj 

STOP  I 

1?  XI  =i-0/x  ; 

BETAl  =1.0/(PM«X>  i 

PSR  = 0.0 

PSRC  = CMPLXfO. 0,0.0)  i 

DC  15  1=1. K2 

MR=K1-I  i 

I«  =?«NR+i  I 

PSR  =1  .0/( IA*XI-PSR)  ; 

PSRC=1 .0/(1 A*BETAr -PSRC) 
if(Nr.gt.l)  go  to  is 

PSTR(NR)  = CMPLXCPSR.  0. 0)  j, 

PSlRr(NR)  sPSRC  . j: 

i5  Continue 


CHIO  = COS(X) 
PSID  = SIN(X) 


DO  50  I =1  ,L 


*»*  CALCULATE  MIF  SERIES  COEFFICIENTS  s** 


I] 


PSR  : REAL(PSIR(  I ) ) H 

II  = I + 1 li 

lr( I .GT.2)  GO  TO  30 

IF(  I .EO.?)  GO  To  20  j' 

pSI  (.1)  =PSR»PS10  II 

PSi  (2)  = REaL(PS1R(2)  )*pSl  (1  ) |1 

CHi (1 )=Xi*CHiO  + PSIO  i 

rMl(?,  ::  3. 0*XUChI  (1  )-CHlO  |- 

9CA(1)=  -XI  ♦ 1.0/(XI-CMPLX(0.0,-1.0))  i 

20  PsI  (.1)=PSl  ( I ) 

CHI (3)=  CHI  ( I ) 

GO  TO  .15 


10  PSI (3)=PSR«PSI(2) 

Ch1<3)  = (2*1-1 )*xI«ChI <2)-CHi(l) 

35  7ETA  5 CMPLXIPSI (3) .CHI (3)) 

ACAPA  = 1.0/PSR  -I*XI 

ACAPB  = l.O/PSIRCd  )-I»BETAI 

RCA (5  ) = -Il*Xl  + 1 .0/1 Il*X I-BCA( 1)  ) 

A(I):  PSI(3)*(  ArAPB-PM*AcAPA  ) / (ZETA*(AcAPB  - PM  «))  CA  ( 1 ) ) ) 
B(I)  = PSi (3)*(ArAPA-PM*ACAPB)  / < ZE TA • ( BC A ( 1 ) -PM» AC  APB  )) 
BCA (t  ) = BCA( 2) 

ICONST  = 2*1  +1 

EXKRO  rICONST  » REAL ( A ( I ) *6 ( I ) ) * EXKRO 

ScKrO  = IcOnST  •(rABS(A(l  )*«2)*CA0S(B(I  )»•?))  ♦ SCKrO 


IF(  I.LT.3)  go  TO  50 
PSl(l)  =PSl(2) 

PSI {?)=  PSI (3) 

CHl ( I ) = CHI (2) 

CHI  ( •’)=  CHiO) 

52  CONTINUE 


EXKRO  = EXKRO  • 
SCKRO  = SCKRO  • 


FIID  = 
LAST  ; 


f2.0/(X**2 

L 


COEPH 

COEPH 

• SCKRO  )) 
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t ‘ »■»«•%<-  s*  ■ ■ * 


ooo  o ooooo 


C 

LMAX  = L • (L+l)/2 

L = length  of  MIE-?  SERIES  ON  A AND  B 

***  CALCULATE  LENGTH  OF  LEGENDRE  SERIES 

IF (X.LE. 15.0)  GO  TO  49 
less  = 21. 0*EXP(0. 040611  • X) 

GO  TO  51 

'^9  LESS  = 2.0  • X * 11.0 

51  IF(LESS.GT.?*L)  LESS  = 2*L 

THIS  IS  the  beginning  of  the  New  recursion  scheme  i 

(I).  COLLINS.  RRA,  7/7/77) 

C i 

MAX  = LESS  ♦ 1 j 

no  52  N=1,MAX 

52  SUMW(N>  = 0.0  I 

DO  100  k=i.last  ! 

Y = K I 

no  100  JsK.last  I 

X 5 J 

hWKJ  - G{K,J)  j 

NsJ-K  ] 

CKJ  = C(K,J)*2.0/M  .0+IFIX(K/J)  ) ] 

IF  (J.EQ.K)  WJKO  = 1.0/(2*Y+i) 

IF  (J.GT.K)  WJKO  = WJKO*X*(  2«X-2*Y-1) /(  (2*X  + 1 )*  (X-Y  ) ) 

WB4  = WJKO 

M = J-K+i  1 

SUMW(M)  = BWKJ*W84*CKJ+SUMW(M ) 

NL  = J-K+3 

NU  r J+K+1  ! 

TOO=2.0/(1.0*IFIX(K/J))  ' ! 

DO  100  M=NL.NU.2  1 

N=M-1  j 

Z = M-1  ] 

WB4  = WB4«(X^Z-Y-l)i(Y-*-2-X-l)*{X  + Y + Z)*(X-fY-Z«-2)/  j 

« ( ( X + Y-Z  + 1 ) *(  X+Y*Z-*^1  )*(  X+Z-Y)*!  Y+Z-X)  ) 

SUMW(H)  = SUMW(M)  + BWKJ»WB4*C(K , J )«T00 
100  CONTINUE 

no  120  Nsl.MAX 
120  SUMV(N)  = 0.0 

no  200  Krl.LAST  1 

Y = K ; 

no  200  JsK.LAST  j 

X = J 

BVKJ  = F(K,J) 

IF  (J.EO.K)  VJKO  = 3*Y*  ( Y*1 )/ (2*Y  + 1 ) Ll 

IF  (J.GT.K)  VJKO  - VJK0*(2*X-2#Y>3)»{X  + 1 )/{ {X-Y)*(2«X*1)  ) 

VB4  ! VJKO 
N : J-K+2 

SUMVIN)  = SUMV(N)+BVKJ*VB4*?.0/(1.0+1FIX{K/J)) 

IF  (K.EO.1)  go  to  200 
NL  = J-K*A 
NU  = J+K 

TOO=2.0/(i.O*IF1X(K/J))  f 


o oo  on  ooo  ono  oon 


!)0  190  MsNL.NU.i? 

Z = M-1 

VB4  : v84*(2*Z  + 1 )*(X  + Z-Y)*'Y  + Z-X)*(X  + Y*Z-*-l  )*(X  + Y-Z  + 1)/ 

S (( 2*Z-3) «<  x*Z-Y-l) *( Y + Z-X-l) «( X+Y+Z) *(  X+Y-Z+Z) ) 

SUMV(M)  = SUMV(M)  + BVKJ«VB4*T00 
lyO  CONTINUE 
?no  continue 

00  ?5o  N=1.MAX 

'^50  COEf(N)  = FUD«(SUMW(N)+SUMV{N)  ) 

c 

c this  is  the  end  or  the  new  recursion  scheme 

ir(PARlS.EO.O.O)  00  TO  402 

95:6  FOLD  LEGENDRE  COEFFICIENTS  WITH  PARTICLE  DISTRIBUTION  »** 

IF(KC.GT.i.AND.KC;.LT.NRS)  go  to  310 
IF (KC .EQ.NRS)  DELR=(RI<KC)-RI(KC-l> )/?.0 
IF(KC  .EO.  1)  DEI.R=  (RI(KC-t-l)-RI(KC))/2.n 
Go  TO  330 

310  DELP=(R1(KC+1)-RI (KC-l) )/2.0 
130  EDRR  =FNR(KC)»nFLR6Rl (KC)*RI(KC) 

CEXT  = CEXT  + EXKR0»EDRR 
SCOE  = SCOE  +SCKR0  * EDRR 


'102  WRlTfc(6,680)  PM  , X YZ  , SCKRO  , L ESS 
WRITE(6,662)  (COEF(N) ,N=1 .LESS) 

PUT  COEFFICIENTS  ON  LIBRARY  (TEMP  ON  NFWLEG,  COPY  TO  LEG) 

WRITF(t2)XYZ,PM.SCKR0.LESS 
WRITEd?)  (CoEFIN)  .N  = 1.LESS) 

IF  (PARIS. EQ. 0,0)  GO  TO  570 


no  450  MMsI.LESS 

COEFXS(MM)sCOeFXS(MM)+COEF(MM)*SCKRO*PI 1#RI <KC)»R1 (KC) 

I :*£NR(KC)«DELR*1  .OF-08 

450  continue  !| 

EnT0T=ENT0T  ♦ ENR(KC)«DELR  • j 

500  continue 

l| 

WRITE(6.7?0)  ENTOT,  LESS 

WRl  TF  (6,555)  (COEFXS(N)  ,N=1.LESS) 

555  FORMAT</, 'MACROSCOPIC  X-S  COEFF  I C I ENTS  ’ . / , ( 1P6E15 . 6 ) ) 'j 

CrXT  r PI  MCEXTbI  .OE-08  i 

SCOE  = PI  I • ScOE  # l.OE-08  i 

A80E=cEXT-SC0E  i 

WRITE(6.710)  CEXT , SCOE. ABoE 


r>  r>n 


5;0  L=p 

560  REAPdl  .END  = 590'  XV7.PM,SCKR0.LESS2 
REALUll)  (COEr(N)  .N  = l .LESS?) 

WRITE(12)XYZ.PM,SCKR0,LESS? 

WRITE (12) (COEF(N) .N=t .LESS?) 

L=L+1 
GO  TO  5>P0 
590  REWrwn  11 
REWIND  1? 

WR1TE(5.595)  L 

595  FORMAT!/. 'READ'  .I?5,  ' SETS  OF  COEFFICIENTS  FROM  LIBRARY  • 

1 'FOR  update. ' ) 

C 

IF{PARTS.E(3,0.0  )STOP 
C 

C OUTPUT  OPTIONS 
C 

599  CALL  OUTOP 
C 

500  FORMAT(6FlO.O, 15) 

^10  FORMAT!/. 2X. * INPUT  D AT A » , // • ?X . • RM I N = F 1 0 2X M ! CRONS ' . 

1/.2X.'RMAX  = ' .F1 0 .2, 2X. 'MICRONS' ./2X. 'PARTS  = '.FlO-?./. 
??X.'PM  = ' ,2E15.6./.?X, 'WAVEl  = ' . F 10 . 2 , ?X , ' M I CRONS ' . 

3/.?X. 'NEO  = ' , 15) 

620  FORMAT!4F10.0) 

630  FORMAT!/. ?X, 'PI  = ' . 1PE15 . 6 . / , 2X . ' P?  = ' . 1PE15 . 6 , / , 2X . 

l'P3  = ' .1PE15.6./2X, 'P4  = '.1PE15.6) 

640  format!/. 2X, 'PI  = ' ,1PE15.6./,2X. 'P2  s »,1PE15.6./. 

12X,'P3  = • ,1PE15.6./,?X, 'P4  = '.1PE15.6) 

650  format!/, 2X, 'NRS  = '.I5./) 

661  FORMAT!/, v3X. 'RADIUS' .6X. 'DISTRIB.  DENS  . ' / ! j P?E15 . 5 ) ) 

662  Format!/, 'LEGENCRE  coefficients,  single  SIZE'/!lP6El5.6)) 

670  FORMAT!/, ?X, 'DIMENSIONS  TOO  SMALL.  L = '.15) 

660  FORMAT!/, 2X, 'PM  = ' ,2F15. 6, 2X  . ' ALPHA  = '.FIS-Z*/. 

IPX.'SCKRO  = ' ,1PE15.6.2X, 'LESS  = '.15./) 

700  format! I10»5X,1P?E15. 7) 

710  FORMAT!/, 2X, 'MACROSCOPIC  EXTINCTION  CROSS  SECTION  = '.1PE12.5. 

1 / ,?X. 'MACROSCOP  IC  scattering  CROSS  SECTION  s '.E12.5. 

2 / .?X. 'MACROSCOPIC  ABSORPTION  CROSS  SECTION  s ',E1?.5.//) 
/20  FORMAT!/. 2X, 'AVERAGE  LEGENDRE  COEFFICIENTS',/. 

1 2X,'ENTf)T  = ' .1PE15.6.2X. 'LESS  = ',15,/) 

STOP 
END 

subroutine  OUTOP 

COMMON/POUT/RMIN.RMAX .PM.WAVEL.NEO.PI ,P?.P3,P4.NRS.LESS. 

1 EN.COEFXS(500) 

DATA  IA0,1A1,IA2.IA3. I A 4/4HXS0R » ^HXSL A , 4HSM00 . 4HAG1 N • AHTRuN/ 
DATA  II . 12. 13. I4/4»0/ 

COMMON  ALPHA 
COMMON/OUT/Il,  12. 13.14 
complex  pm 

other  output 

REAn!5.302) lOUTl, I0UT2, I0UT3. I0UT4 
302  F0RMAT(4(tx.A4)  ) 

C 
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CHECK  HOLLERITH  VALUES 


WRITE(6.303) lOUTl , I0UT2. IOUT3. IOUT4, I AO . I A1 , I A2 . I A3 . I A4 
303  FORMATf/, 'OPTIONS  CHOSEN:  ’ , 4 ( 1 X . A4 ) , / , • A VA  IL ABLE  OPTIONS:'. 

1 5(1X.A4)) 

cccccccccccccccccccccccc 


11=1  IF  x-S  SET  FOR  ANISN  (OR  OTHER  OAK  RIDGE  CODE)  PRODUCED 
IF( lOUTl.EO. 1 AO .OR. 1 0UT2 . EO • I AO . 0? . I0UT3,EQ. I AO. 

1 OR. I0UT4.E0. lAO)  Il=l 
I1WAS=I1 

C tl=2  IF  X-S  SET  FOR  DTF  (OR  OTHER  LOS  ALAMOS  CODE)  PRODUCED 
IF ( lOUTl .EQ. lAl.OR . I0UT2.E0. 1 Al . OR . I 0UT3 ,EQ . 1 A1 . 

1 OR. I0UT4.E0. I Al ) Tlr2 
C 12=1  IF  SMOOTHING  TpsTED 

IF( 10UTl.EQ.IA2.OR. 10UT2.EQ. IA2.0R. I0UT3.E0. IA2. 

1 OR. I0UT4.EQ.IA?)  12=1 
C 13=1  IF  REFITTING  OF  COEFS  TO  F TESTED 

IFdouTl.EQ.  1A3.0R.  1OUT2.E0.  I A3 . OR  . I 0UT3  .EQ . 1A3. 

1 OR . I0UT4 .EQ. IA3)  I3=l 
C 14=1  IF  truncation  tested 

IF( IOUTI.EQ. I A4.0R. I OuT 2 . EQ . I A4 . OR . I 0UT3 . EQ . I A4 . 

1 OR. IOUT4.EQ. 1A41  14=1 

C 

IF(  Il.EQ.O)  GO  TO  288 
IF( IlMAS.EO.l.AND.ll.EQ.2)  GO  TO  2^0 
GO  TO  288 
?;o  Ii=i 

WRITE(6,272) 

272  format:/, 'CANnO'  PRODUCE  O-R.  AnD  L.A.  FORMAT  XS  ' 

1 'SETS  IN  THE  SAME  RUN.  O-R.  (E.G..  ANISN)  CHOSEN.') 

287  FORMATC  output  I ND  I CFS  : ' , I 5 , 5X  . 3 1 5 ) 

280  FORMAT! • OUTPUT  INDICES:  '.SX.AIS) 

288  IF( II .EQ.O)  go  to  305 
IF( II .EQ. 1)  GO  TO  285 
Do  280  MM=1,LESS 

C0EFXS(MM)  = C0EFXS(MM)/(2.''FL0AT(MM)-'1.  ) 

280  CONTINUE 

WRI TE(6 ,289)1 1, 12,13, 14 
GO  TO  ?81 

285  WR1TE(6,287)I1, 12,13, 14 

281  7ERO=0.0 
M0RDER=2.*ALPHA+3. 

WRITE(15,4) 

4 FORMAT'/, 'MAC  DATA  PRODUCED  BY  MIELEG') 

WRITE(15,291)R'MIN,RMAX,PM,  WAVEL,NE0,P1  . P2 , P3  iPA,  NRS , LESS  . EN , MORDER 

291  FORMAT! /, 'PARAMETERS' ,/ ,5E15. 6. /, 14 ,4E15. 6, / . 21 5, E15. 6, 15) 

WRi TE( 15, 292) zero, ZERO, COEE XS ( 1 ),C0EFXS(1  ), 

1 : ZERO. ZERO. ZERO. COEFXS! MM) .MM=?,LESS) 

292  FORMAT!/, '6E12. 5 FORMAT  CROSS-SECT  1 ONS ',/,' 14«* ', /2F 12 . 1 . 

1 ?ei?.5,/, (3F12.1.E12.5)  ) 

IF! II .EQ. l)WRI IE!6,293) 

293  FORMAT !/. 'MACR  X-S  DATA  WRITTEN  ON  MAC  FILE  (15).-  OPTION  XSOR * ) 
iFdl.Fo.?)  WRITB(6,294) 

^94  FORMAT!/, 'MACR  X-S  DATA  WRITTEN  ON  MAC  FILE!15).  - OPTION  XSLA') 

C 

T05  CALL  TESTAN!0.rOEFxS.LESS) 

RfTuRN 

END 
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subroutine  TESTANI I , AN.LESS) 

DIMENSION  X(l8l) ,AN(500) ,AN?(500) , AN3 ( 500 ) , P ( 500 ) , 

1 PP(500).p3(181).F(181).F2(18l),r3(18l) ,N(iel).TW(101). 

2 FFlTdSl) 

COMMON  ALPHA.Pl , IP. 1P2. IP3 
COMMON/OUT/Il,I2.l3.14 

if(lrss.lt.5oi)  go  to  10 
WRITE(<»,5)  less 

5 FORMAT! /. 'LESS  = *,14.'  MUST  CHANGE  DIMENSION  IN  TESTAN  aND  LEGCHK 
1’  ) 

STOP 

10  IM  1 .GT  .1  ) GO  TO  4? 

TlsPl/iao. 

DO  30  J=l,181 
T=Tl*(J-l) 

X(J)=C0S(T) 

30  CONTINUE 

PRODUCE  leg.  POLyS  AND  STORE  WITH  ALL  ORDERS 
FOR  ONE  ANGLE  ON  RECORD 

1P  = 8 
IP?  = 9 
REWIND  IP 
REWIND  IP2 
DO  41  J=l,181 
P(1  ) = l.O 
P(2)=X( J) 

DO  40  L=3.LESS 
FL=L 

P(L)={ (2.0«FL-3.0)#(X( J)«P(L-1) )-(FL-2.0) 

1 •P(L-2) )/(FL-1.0) 

40  CONTINUE 

WRITE! IP) (P(L) .L=1 .LESS) 

41  CONTINUE 
REWIND  IP 

PRODUCE  leg.  POLYS.  AND  STORE  WITH  ALL  ANGLES  FOR 
ONE  ORDER  ON  RECORD 

lP3rl0 
REWIND  IP3 
DO  135  J=1.181 
:.35  P(J)«1  .0 

WRlTt;(  1P3)  (P(J)  .J  = 1.181) 

00  l36  J=l.l8l 
P2(J)  = P(J) 

1T6  P(J)  = X(J) 

HRITE(IP3)(P(J).J=1.181) 

DO  I4i  L=3.LESS 
FL  = L 

DO  140  J=l«l8l 
P3( J)=P2(J) 

P2( J)=P< J) 

P(J)  = { (2.»FL-3. )*(X( J)*P2( J) )-(FL-2. )*P3(J) )/(Fl-1.  ) 

1 40  CONTINUE 

WRITE(IP3)(P( J) .J=l.l8i) 
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141  continue 

REWIND  IP3 
GO  TO  46 


42  IF(LESS.LE.LESWAS)  GO  TO  46 


L1=LESWAS+1 

no  45  J=l,181 

REAnnP  ) (P(L)  .L=t  .1  ESWAS) 

00  44  L=L1.LESS 
FL  = L 

P{L)=( (2.0*FL-3.0)«(y( J)*P(L-1 1 >-<Fl-?.0) 

1 *P(1.-2))/(FL-1.0) 

44  continue 

WRITE{IP2)(P(L) .L si, LESS) 

45  CONTINUE 
I tEmPsIp 
1 P= IP2 
IP?=ITEMP 
REWIND  IP 
REWIND  IP2 

! WMsl.ESWAS-2 
no  144  L=1,LWM 
1 44  READ( IP3) 

READ( 1P3) (P2( J) . Jsl.lSl) 

READ(IP3)(P<J),Js1,181) 
no  14q  L=L1,LESS 
FI.sL 

DO  146  J=1,181 
P3{ J) sP2( J) 

P?( J)sP(J) 

P( J)  = ( (2.*FL-3, )o(X( J)*P2( J) )-(Fl-2. ) *P3 ( J ) ) /( FL-1  . ) 
1 46  CONtI  nUE 

WRITE(1P3)(P(J) »J=l.l8l) 

148  CONTINUE 
rewind  IP3 

46  leswas  = less 

CAUL  CALcF(AN.LESS,F) 

IFd.QT.O)  RETURN 

evaluate  f using  truncated  series 

IF{14.EQ.1)CALL  LEGcHKI AN.LESS) 

JSE  smoothed  values  OF  F TO  GET  LOwER'ORDER  FIT 

IF< T2.EO.O)GO  to  171 

CALL  SMOOTH{F,N,TW,FFIT,NFl T) 

L=1 

163  SUMS  0.0 

REAn( IP3) (P{J), J=l, l8t) 
no  165  J2sl  ,NF1  T 
NPT  =N(J2) 

X J=X(  NPT) 

S IN JsSORT (1 .-XJ«X J) 

SUMsSuM>SINJ»FF  IT  {J2)  -TWI  J2  )«P(  NPT) 

165  CONTINUE 
EIMtL-1 

SUMsSuM*{2.*ELM>1 .)*2.*PI 
AN2(I  )sSUM 


oonoo  non  oon  ooono 


l =L+1 

IFd.LE  .LESS)  GO  TO  163 
REWIND  IP3 

WR1TE(6,170)(AN?(L) ,L=1.LESS> 

170  FORMAT! /. ’COEFS  FROM  SMOOTHED  F • , / . ( 1P6E1 5 . 6 ) ) 
CALL  CALCF( AN?,LFSS,F?) 

IF(M.EQ.1)CALL  L EGCHK  ( AN2  - LESS  ) 


INTEGRATE  F OVER  ALL  ANGLES.  RECALCUl ATE  AN(L)  VALUES 

(F  not  SMOOTHED) 


1/1  1F( 13.EQ.0)RETURN 
L=1 

WR1TE(6,62) 

62  FORMAT!/. 'RECALCUL ATED  COEFFICIENTS  - OPTION  AGIN.'./. 

1 '!SlMlLAR  TO  1 npG  SMOOTHING.)') 

63  SUMrO.O 

REAn! IP3) (P( J) . J=l. l8l) 

DO  66  J=l,180 
COSJ= !X(J)+X! J+1) )/2. 

SINJ=S0RT!1.-C0SJ*C0SJ) 

Fj=!F(J)+F! J+1) )/2. 

ELM  = L-1 

PWT  = < 2 - »ElM+1  . ) «(P(  j)  +P  ( J+1 ))  /2  . 

SUM  = SUM'-S1NJ*F  J*TlePWT 

6^  continue 

SUM=  SUM*2.*PI 

WRl TE (6  ,7o)L. SUM, AN!L  ) 

/O  FORMAT! 'L. SUM, AN!L)  = ' . 1 5 , 1P2E1 5 . 6 ) 

AN3(L )=SUM 
IFd.NE-O)  RETURN 
L=L+1 

1F(L .LE .LESS)  GO  TO  63 

REWlNn  1P3 

call  CAI.CF!  AN3.LESS) 

TRY  TRUNrATING  THE  RECALCULATED  SERIES 

IF( 14  ,EQ. 1) CALL  LEGCHK! AN3.LESS ) 

RETURN 

END 

SUBROUTINE  LEGCHK (AN. less  ) 

COMMON  ALPHA, PI, IP, 1P2.IP3 

DIMENSION  FLO!101 ) .LORDER!/ ) , AN!5  00  ) .LTEST! 6) ,P(500  ) 

CHECK  FOR  LOWER  LEGENDRE  ORDER  ACCURACY 
WRITE!6,2)LESS 

2 FORMAT!/. 'TRUNCATE  TO  LESS  THAN', 14,'  TERMS.  !OPT10N  TRUN) 
1 'TEST  FOR' ) 

LORDER! 1)=0 

CCCCCCCCCCCCCC 

FIND  TRUNCATION  TEST  POINTS 


o o n o o o 


no  5 1=?,6 

5 LTESK  n = l 
ITEST=1 

LTEST(1)=2.*AlPHA43. 

WRITE (6,8)L1EST(1  ) 

8 EORMAtMS,’  {?«ALPHA  + ?)  • ) 

C 

c look  for  local  minimum 
c 

L?=LEsS/? 

6 no  10  L=L2,LESS 

If  ( AN(L ) .GT. AN(L-1 > .ANO. ANf L) .GT. AN(L-2)  ) GO  TO  1? 

10  continue 

Go  TO  15  . h 

1?  I TESTrl TEST+1 

LTEST(ItesT)=l 

WRITE(6.13)L  . , 

l3  FoRMAT(I5,'  (LOcAL  MINIMUM)') 

IF( ITEST.E0.5)  GO  TO  15 
L?=L+1 
GO  TO  ft 

j 

LOOK  FOR  first  COEFF I C I ENT . L T . AN < 1 ) j 

15  no  20  L=2.LESS  i 

IF(AN{1  ).LT.AN(1)  ) GO  TO  22  I 

^0  CONTINUE  I 

?2  ITEST=ITEST+1  i 

LTESTI ITEST)=L 

WRITE(6,23)L  S 

^3  F0RMAT(l5.’  (LT. FIRST  COED')  | 

PUT  TEST  POINTS  IN  ORDER 

1F(  ITEST.EO.l  ) GO  TO  27 
ILOOP=0 

I TESTM= ITEST-1 
2A  IF( ILOOP.GT .15)  SIOP 
DO  26  1=1 , ITESTM 

1F(lTEST( I ) .lE.lTESTI  I + l)  ) GO  TO  26 
I TeMP=LTeST( I ) 

LTEST ( I )=LTEST( 1+1 ) 

LTEST( 1+1 )=ITEMP 

iloop=iloop+i 

GO  TO  24  . 

^6  CONTINUE  I 

C i 

C calculate  resulting  phase  FUNCTIONS  j 

^ ] 

27  00  2A  I=1,ITEST 

LORDER  ' I + 1)=LTEST(  I ) 

^8  CONTINUE  i 

C I 

no  125  j=1.101 
125  FLO(J)=0.0 


C 


L?=L0R[JER(L0  + 1) 

Do  155  J=l.l8l 

RtAD(  IP)  (P(L)  »l.  = l .LESS) 

no  l5n  l=Ll,L2 

n o(j)=rLO(j)^AN(i  )»p(L)/(4.«pi ) 

’50  ("ONTInUe 
CONTINUE 
REWINr  IP 

WR1TE(6,156)L?.  (F|.0(J).  J=l.  181) 

156  FORMAT!/, ‘L  = ' . M , / , 1PE15 . 6 . / . ( t P6E1 5 . 6 ) ) 

1.10  continue 

RETURN 

ENP 

subroutine  CALCF( AN.LESS.F) 

COMMON  XYZ.PI .IP, IP?. IP3 
DIMENSION  P(500).F(l8l) ,An(500) 

DO  55  J=1  ,181 
F(J)=0.0 

REAP!  IP)(P(L)  .Lrl  ,LESS) 

PC  50  L=i,UESS 
F (J  ) = F(  J)+AN(L)  *P  (I-  ) 

50  CONTINUE 

F(J):F(J)/(4.*Pn 
57,  Continue 

REWINO  IP 

WRITE(6.60)(F{J),J=1,i8i) 

60  FORMAT!/, 'F  FROM  LEO.  POLYS.,  0 TO  180  DEGREES *./, IPE 15 . 6 , 
1/,!1P6E15.6)) 

RETURN 

End 

subroutine  SM00TH!F ,N.TW,FFIT,NFIT ) 
dimension  F!l8l ) .N(iei) ,TW!181) ,FfIT!1B1) 

T1=3. ’415926536/180. 

WRITE!6.245) 

Z'lS  format  !/, 'SMOOIHING  PARAMETERS’) 

REAn!5.246)THRESHF,WlnTH 
?46  FORMAT!2E10.0) 

WRi TE !6 ,247 )ThReShF ,W iDtH 
?47  FORMAT!/, • THRESHF . W I DTH= ’ ,1P2E15,6) 

DO  ?48  J=2.90 
F1=F(’ )*THRESHF 
IF!F! J) .LT.Fl)  GO  TO  ?49 
.’T8  continue 

WrItE!6,250) 

250  FORMAT!/, 'NO  VALUE  OF  PHASE  FUNCT  IN  FIRST  90  DEG  LESS  THAN’,/, 
1 ’F(’  )»THRESHF  to  begin  SMOOTHING') 

GO  TO  61 
2a9  JwlDErWIDTH 
w1()Th=JWI  OE 
NlNT=l8n. /WIDTH 
PO  251  Irl.NiNT 
JLASTrUjMlDE+l 
IF! J.LT.JI AST)  GO  TO  255 

251  Continue 
WRITE!6,252)  J 

252  format!/. ’THRESHOLD  ANGLE’. 14, • NOT  FOUND’) 

GO  TO  61 


2bb  NF! T= J+I+NINT-I 
WRITE(6.253) 

Pb-i  format(/,'rep.  pt.  index  width  F*) 

N(1  )rl 
TW(1)=T1/2. 

FF1T(1)=F(1) 

DO  256  K=2.J 

N(K)=K 

TW<K)=T1 

256  FFIT{K)=F(K) 

WRITE(6,254)(N<K) ,TW(K).FFIT(K).K=1,J) 
FORMAK  I5,2E15.6) 

N(  ) = ( JLAST  + J)/2 

FFI TJ=0.0 
JM  = Jl  aST-1 

IF(N(  J-1.)  .EO.  J)N(  J-t-DrJ  + l 
IF( J.EO. JM>GO  TO  357 
DO  257  J2=J+1»JM 

257  FFlTJrFFITJ+F(J2) 

<57  FF1TJ=FF1TJ+F( JLAST>/2. 

P01NTS  = FL0AT( JLAST-J)-0 .5 
FFIT(J+1)=FFITJ/P0INTS 
TW( J+1 )=POINTS*Tl 

WRITE(6,?54)N( J + 1 ) ,Tw( J+1) .FF1T( J + 1 ) 

JMln=jLAST-JWlDE/2 

DO  260  KsJ-*-2,NFIT 

JMlDsjMID-fJWlDE 

N(K)=JMID 

jfirst=Jlast 

JLAST  = JF  IRST4.JWI  nE 

FFITKr  (F( JFiRSr )+F( JLAST) )/2. 

TFIJWIDE.EO.I)  Go  TO  259 

JP=JFIRST+1 

JMsJL AST-1 

DO  25fl  JsjP.JM 

FFI  TK=FFI  tK+F  ( J ) 

250  CONTINUE 

''59  FFIT{K)=FFITK/FL0AT(JWIDE) 
TW(K)=T1»FL0AT( JWIDE) 

WRI TE(6 .254 )N(K  ) , TW(K  ) , FF 1T(K  ) 

’60  CONTINUE 
51  RETURN 
END 


APPENDIX  B.  LEGEDIT  UTILIZATION  INSTRUCTIONS 


The  LEGEDIT  code  lists,  orders,  and  deletes  sets  of  microscopic 
Legendre  coefficients  on  the  LEG  file.  It  can  be  Included  in  a run 
after  MIELEG  to  order  the  file  produced  from  a pre-existing  file  and 
the  new  MIELEG  contributions,  or  It  can  be  run  by  Itself. 


B.l  LEGEDIT  Input  Description 


MODE  (15) 

LEGEDIT  Input  data  consists  of  a single  Integer  variable, 

MODE.  The  options  are: 

MODE  - 0 Identify  entries  on  files  LEG,  NEWLEG,  BACKLEG,  OTHERLEG 
LEG  - Normal  Legendre  coefficient  file 
NEWLEG  - Temporary  production  LEG  file  during  MIELEG 
run. 

BACKLEG  - Back-up  LEG  file,  produced  at  beginning  of 
MIELEG  run. 

OTHERLEG  - Overflow  file  for  other  refractive  Indices,  3 
Contains  entries  edited  from  LEG  file  which 
did  not  have  same  refractive  Index  as  ori- 
ginally found  on  file. 

^K)OE  - 1 Ordering  mode.  Puts  entries  In  order  of  ascending  size 
parameter,  a.  Deletes  old  entry  If  value  of  a within 
0.01  of  new  entry.  Writes  on  OTHERLEG  any  entries  with 
refractive  index,  3,  different  from  that  originally  on 
file. 


MODE  - 2 
MODE  - 12 
MODE  - 21 
MODE  - 212 


List  mode.  List  LEG  file 
Order,  then  list . 

List,  then  order 

List,  order,  then  list  again. 


Deletes  selected  entries. 


67 


laiMiKiiteMiiyii 


MODE  - 3 


ALFOUT  (ElO.O)  Ignored  unless  MODE  - 3 

Up  to  10  selected  deletions,  input  one  per 
record.  Entries  will  be  deleted  If  size  parameter, 
a,  within  0.01  of  ALFOUT. 

B.2  LEGEDIT  Output  Description  - Printout 

Identification  mode  0 prints  size  parameter  (a, ALPHA) , complex 
refractive  Index  (3, PM)  cross-section  ratio  (K,  SCRRO),  and  order 
(N,LESS)  for  each  entry  of  LEG,  NEWLEG,  BACKLEG,  OTHEKLEG.  Ordering 
mode  1 (or  the  ordering  part  of  12,  21,  or  212)  lists  the  size  parameter, 
complex  refractive  Index,  cross-section  ratio,  and  Legendre  order  of 
each  entry  on  file  LEG.  This  list  Is  produced  once  before  the  size 
parameter  ordering  and  deletions  occur,  and  again  after  they  occur. 

List  mode  2 (or  the  listing  part  of  12,  21,  or  212)  lists  the 
size  parameter,  refractive  Index,  cross-section  ratio,  and  Legendre 
order  of  each  entry  on  LEG,  followed  by  the  set  of  coefficients  asso- 
ciated with  these  parameters.  Mode  3 (selected  delations)  operation 
is  followed  by  a list  of  the  amended  file. 

B.3  LEGEDIT  Output  Description  - Files 

The  purpose  of  the  ordering  mode  of  LEGEDIT  is  to  produce  a 
LEG  file  In  ascending  order  of  the  size  parameter.  In  addition,  the 
LEGEDIT  Job  (OR  MIELEG  followed  by  LEGEDIT)  saves  a backup  file, 

BACKLEG,  which  Is  a copy  of  the  original  LEG  file.  The  updated  LEG 
file  Is  produced  on  a file  named  NEVRiEG,  so  an  Incomplete  Job  might 
have  the  required  data  on  NEWLEG.  A successfully  completed  LEGEDIT 
Job  would  produce  Identical  LEG  and  NEWLEG  files. 


B.4 


LEGEDIT  FORTRAN  Listing 


nooonnnnn  oonn 


lEGEDIT  --  code  TO  FPIT  I IBRARY  OF  LEGENDRE  COEFFICIENTS  FROM  MIELEG 


HEAr(5»10)  mode 
10  FORMAT! 15) 

MODF  determines  EDITING  TYPE 

=0  To  READ  FILES  11,12,13.14  ( LEG . NEMLEG ,8 ACKLEG ,0 THERLFG ) 

-1  TO  order  and  delete  DUPIICATES  (ALWAYS  DELETE  WITH  ORDER), 

WRITE  entries  with  DIFFERENT  REFR.  IND.  ON  OTHERLEG  ,1 

= ? TO  LIST  EN1  IRE  FILE  ij 

=12  TO  ORDER  THEN  LIST 
=21  TO  list  then  order 
= 212  TO  list,  order,  THEN  I 1ST  AGAIN 
C = 3 TO  delete  selected  ENTRIES 

C 

REWIND  11 
REWIND  52 

IF(M0DE .E0.3)CALL  SCRUB 

if(mode.eo.o)  call  readit 

IF(M0rE -EQ.D  call  ORDER 
IF(M0DE.EQ.2)  call  LlST(ll,12) 

IF(M0DE.NE.12)  GO  TO  21 
CALL  ORDER 
call  LIST(12,12) 

?1  IF(M0DE  .NE.21 ) GO  TO  ?1.? 

CALL  LlST(ll.ll) 

Call  oRder 

21?  1F(M0DE  .NE.212)  STOP 
CALL  LlSTdl.ll) 

CALL  ORDER 
call  LIST(12.12) 

STOP 
fnd 

subroutine  readit 

C tile  11  (LEG)  - NORMAL  LIBRARY  (SINGLE  IND.  OF  REFR.) 

C file  12  (NEwLEG)  - NEW  ADDITIONS  TO  LIBRARY,  MAY  NOT  BE 
C OPEN  (IT  LAST  addition  EXITED  NORMALLY) 

c file  i3  (BAckleg)  - leg  backup 

C fILE  14  (OTHERLEG)  - RECEPTACLE  FOR  DATA  WITH  OTHER  INDICES  OF  REFR. 
complex  pm 
DO  50  1=1.4 
IN  = 10  ♦ I 
WRITE(6,19) IN 

19  FORMAT!/, 'CONTENTS  OF  FILE'-I3) 

•'0  READ!  IN,  END  = 50  ) ALPHA,  PM,  SCKRO,  LESS 
WR I Tp( 6.15) ALPHA , PM, SCKRO, less 
15  F0RMAT(/,F7.2.5x.2F7.2.5X,1PE15.6, 110) 

REA  T( IN) 

Go  TO  20 

50  continue  li 

REWIND  IN  il 

RETURN  i 

End  I 

subroutine  order  I 

complex  PM(  200) ,P  .PlAST  i 

dimension  ALPHA(?nO), less (200  ), lORDER (200  ) , COEF (50D  ) , I TAG (200  ) ^ 


1 .SCKRO(20O) 


REWlNp  14 
I =0 

1 other  = 0 

INEMaS'O 
10  I = I 1 

READ(  11  .ENDsSO)  ALPHAd  ),PH(1  ).snKRO(  I)  .LESS(  n 
READ(  1l ) 

ir(  I.LT  .200)00  TO  10 
WRI TE(6,15> 

’3  FORMAT( /, 'NUHBER  OF  FILES  EXCEEDS  200  . DIMENSIONS  TOO  SMALL 

return 


•jO  ILAST  =1-1 
PLAST=PM( ILAST) 

WRITE (6,55) ( 1 , ALPHA ( I ),PM( I ) .SCKRO( I ) , LESS! I ) . 1 *1 . 1 LAST  ) 
55  FORMAT{/, IIX, 'ALPHA ». 4X, » IND.  OF  REFR . • ,8 X, * X-S  RATIO'. 4X 

1 .'ORofR  ON  LIB  FILE'/ 

2 {14,5X.0pF7.2,5X,2F7.2,5X.1PE15.6,  18)) 


no  57  1=1, ILAST 
I TAG(  T ) =1 
lORDERf  1 ) = I 
57  continue 


ISAVF  = IlAST 
ILASTM= ILAST-1 
58  11  : 1 

‘’R  nO  60  1=11  , ilastm 

IF(PM(1 ).eO.PlAST)GO  to  62 

1F(lTAG(I).EO.O)GO  TO  60 

WRI TF (14) alpha ( 1) .PM( p ,SCKR0(1 ) ,lESS( I ) 

WRITE(14){C0EF(N) ,N  = 1,LESS(  I ) ) 

WRITE(5,61)I 

61  FORMaTC  entry  no. ',14,'  WRITTEN  TO  FILE  OTHERLEG  SINCE  IND 
1. ' refr  different . ' ) 
iother=i 

I TAG( I >=0 
GO  TO  63 

t2  1F(  APS*  ALPHA!  I )-ALPHA(  Ul) ) .LT.O.ODGO  TO  68 
IF  ( alpha  ( 1) -ALPHA  ( I -t-i)  ) 60. 68, 65 
68  IF( 1TAG(1+1).EO.O)GO  TO  60 
I TAG(  l + DsO 


63  ISAVE=ISAVE-1 

60  continue 

IFdi.GT.l)  GO  TO  58 
GO  TO  70 

65  I TEMP  = IORDER(  I) 
ITAOTsITAGI I) 
ATEMP:ALPHA( I ) 
lORDERI  1 )sIORDER(  PI) 
ITAG(  I ) = ITAG(  PI) 
ALPHAd  )sALPHA(  Pi  ) 
lORDER(Pl)  = ITFMP 
ITA6(  P1)«1TAGT 
ALPHA(  I 'DsATEMP 


IFC 1 .EQ. ILASTM)G0  TO  58 
^6  11  = 1 + I 

GO  TO  59 

70  REWIND  It 

rewind  1? 
rewind  t4 

TNEW  =1 
WhItE(6,71) 

71  roRMAT(//, 'REORDERED  AS  FOLLOWS.  (DELETIONS  SKIPPED.)') 
(•Q  DO  80  I=1,IlAsT 

IF( 10RDER{ INEW) .NE. I ) GO  TO  75 

IF( ITaG(InEW) .EO.O)  GO  TO  73 

READdDA.P.SK.L 

READ! 1 1 ) (COeF (N ) .N=1.LESS ( I ) ) 

WRITE (6 .72) INEW.A .P.SK.L 
/?  FORMA T< IA.5x.r7.2.5x»2F7.2-5x,lPEl5.6. 18) 

WR1TE(12)  a, P.SK.L 

WRITE (12) (COEF(N) .N=l ,LESS(  I)  ) 

INEW  = iNEW  + 1 
GO  TO  74 
'3  READ  ( 11 ) 

READ  (11) 

INEW  r InEW  1 
74  IF(INEW.GT.ILAST)  GO  TO  90 
r,0  TO  80 
■’5  READdl) 

REAiXll) 

60  CONTINUE 
rewind  ll 

IF(  INEW.GT.  ILASDGO  TO  90 
IF( INEW.EO. INEWAS)  GO  TO  82 
!NEWAS=INEW 
GO  TO  69 

P.P  WRI  TE(6,84)  INEW 

H4  FORMAT!/. 'CANNOT  FIND  NO.'.U.'  STOP.’) 

return 

OC  J=0 

IF(  IOTHER.EQ.OGO  to  110 
WRI  TE  (6  .92) 

9?  Format!/, 'Entries  on  otherleg') 

95  READ(14.END=110)A. P.SK.L 
JrJ-t.1 

REA D( 14 )( COEF (n).N=1, LESS (1  )) 

WRITE(6,72)J, A. P.SK.L 
GO  TO  95 
1 lO  REWIND  It 
rewind  12 
rewind  14 
RETURN 
End 

subroutine  LISTdN.IOUT) 
coMpi  px  Pm 
dimension  COEF(5oO) 

WRI TE (6,10) 
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i 

1 

i 


j 

I 

I 


IV 

L2. 


f 


10  FORMaK/,  •contents  of  leg  library  file*) 

I =0 

RFAD(  IN.END*30)ALPHA.PM.SCKR0.LESS 
READ(  INXCOeE  (N).N=1.lESS) 

1 = I+1 

IF( IN.EQ. lOuT)  GO  TO  14 

WRl TE( lOUT) ALPHA.PM.SCKRO.LESS 

WRITE(IOUT) (COEF(N) .Nrl.LESS) 

14  WR1TE(6,i5)  1 .alpha, pm, SCKRO, less 

FoRMAT(/.  'Entry  •,  13, 5x.  •ALPHA  = 'F7,2,5x.  'PMs*  .2r7.2,5X, 
1 'SCKROs*  ,iPE15.6.5X,  •LESS=M4) 

WRItE(6,20) {COEF(N) ,N=l,LESS) 

20  FORMAT!/, (1P6EJ5. 6) ) 

Oo  TO  12 

30  rewind  In 
Rewind  iojt 
return 
End 

subroutine  scrub 

DIMENSION  AlFOUT(10),cOEF(500) 
complex  pm 
DO  10  1=1,10 

READ(5,5,END=15)ALF0UT(  I) 

0 Continue 
5 formatieio. 0) 

15  IoUT=l-l 

WRl TE(6 ,1 6) (ALFOUT! 1 ) ,1 =1 ,IOUT) 

•6  format!/. 'SCRUB  LEG  ENTRIES  WITH  ALPHA: '/ !FlO  .2)  ) 

20  READ!  11  ,END=50)  alpha,  PM,  SCKRO,  LESS 
READ! ll)!C0Er!N) .N=1,LESS) 

DO  25  I =1,1  OUT 
DALFsABS!  ALFOUT  ! I )-ALPHA) 

IFIDALF .LT. 0. 01 ) GO  TO  20 
25  continue 

WR I TE!i2) alpha , PM, SCKRO, less 
write !12) !COEf!N) ,N=1,LESS) 

GO  TO  20 
50  rewind  t2 

CALL  LIST112,12) 

STOP 

END 


APPENDIX  C.  MACRO  UTILIZATION  INSTRUCTIONS 


The  MACRO  code  produces  a set  of  macroscopic  Legendre  coeffi- 
cients from  the  microscopic  sets  on  the  LEG  file.  MACRO  Is  essentially 
identical  to  parts  of  the  MIELEG  code,  which  produces  macroscopic 
coefficients  directly  from  basic  parameters  of  wavelength,  refractive 
Index,  and  particle  size  distribution  and  range. 

C.l  MACRO  Input  Description 

TITLE  (20AA)  Any  descriptive  alphanumeric  Information. 

RMIN,  RMAX,  PM(Complex),  WAVEL,  NEQ  (5F10.0,  15) 

RMIN  - Minimum  particle  radius  considered,  microns 

RMAX  - Maximum  particle  radius  considered,  microns 

PM  - Complex  refractive  Index,  3,  composed  of  two 
parts i 

a)  M,  real  Index  of  refraction  relative  to 
the  surrounding  material,  wd 

b)  K,  extinction  coefficient  of  the  particle 
material. 

NAVEL  - Wavelength  of  the  radiation,  microns 
NEQ  - Number  of  equation  describing  particle  size 
distribution,  0 or  1 
0:  N(R)  - P1*R**P2*EXP(-P3*R**P4) 

1:  N(R)  - PI  for  R £ P2 

N(R)  - P3*R**P4  for  R > P2 

PI,  P2,  P3,  P4  (4F10.0) 

Parameters  described  by  choice  of  NEQ,  on  previous 
Input  record. 

Output  Options  4 (IX,  A4) 

Possible  options: 

XSOR  or  XSLA,  TRUN,  AGIN,  SMOO 
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XSOR  - Produces  macroscopic  Legendre  coefficients 

appropriately  normalized  for  use  In  Oak  Ridge 
code  (ANISN,  TDA,  DOT,  MORSE). 

XSLA  - Produced  macroscopic  Legendre  coefficients 
appropriately  normalized  for  use  In  Los 
Alamos  code  (DTP,  TRAMZIT,  TNOTRAN) . 

TRUN  - Truncation  option  to  test  expansion  of  Legendre 

series  to  fewer  than  total  nund>er  of  coefficients 
available . 

AGIN  - Refitting  option  to  reproduce  coefficients  from 
expansion,  testing  consistency  of  method. 

SMOO  - Smoothing  option  to  refit  from  smoothed  data. 

THRESHF,  WIDTH  (2E10.0)  - Smoothing  parameters,  used 
with  SMOO  option. 

See  MIELEG  Utilization  Instructions  for  a more  complete  discussion 

of  output  options. 

C.2  MACRO  Output  Description  - Printout 

The  MACRO  code  prints  out  all  Input  data.  It  then  prints  the 
size  parameters  of  all  data  from  the  LEG  file  appropriate  to  the  wave- 
length and  size  range  of  the  problem,  along  with  the  particle  size,  range, 
and  relative  density  associated  with  each  of  the  data  points.  It  then 
prints  the  complete  LEG  file  entry  for  each  of  these  'data  points. 

If  the  TRUN  option  Is  used,  the  expansion  Is  printed  for  each 
truncation.  If  the  S>KX)  and  AGIN  options  are  used,  the  resulting 
coefficients  are  printed  along  with  the  full  expansion.  The  TRUN  option, 
along  with  S^KX)  or  AGIN  also  causes  the  truncated  expansions  for  these 
cases  to  be  printed. 


,3  . L 


r>  n o r>  n 


LEGENDRE  COEF  ON  FILE  11  (LEG)  AND  PUT  X-S  ON  FILE  15  (MAC) 


complex  P.PM 

DIMENSION  COEF(500  ) , COEFXS  ( 5(10  ) . R I ( POO  ) , ENR  ( 200  ) , ALPHA  ( 200  ) . 
1 SCKRO(200) .SIG(200  ) .LI ( 200 ) . DELR ( 200 ) , T 1 TLE ( 20  ) 

Data  iao.iai,ia?.ia3. taa/ahxsor.ahxsla. ahsmoo.ahagin. ahtrun/ 

DATA  11,12, 13, I4/A*0/ 

COMMON/OUT/Il, 12, 13, 14 
COMMON  ALPHER.PI 
PI =3. iAi5926536 
READ(5.?)  (TITLEd  ),  1 = 1,20) 

2 format(?oa4) 

WRITE(6. 4) (T1TLE( T ), 1=1,20) 

1 FORMAT(/.20A4) 

DO  10  MH=t,500 
■0  GOEFXS(MM)=0.0 


'DETERMINE  DATA  REQUIRED 


RFA[)(  5,  15  )RMI  N,  RMAX  ,PM,  WAvEL,  NEQ 

15  FORMAT(5fio.O,I5) 

WIRI  TE  (6  ,16)RM  IN  ,RmAx.  PM  .WAVEL 

16  FORMAT!/, •RMIN,RMAX='2F7.2» ' MICRONS  (E-04  CM)',/ 

1 'Index  of  refr  = ' 2F7. 2, ' (complex)',/ 

2 'WAVELENGTH  ='F7.2) 

W AVEN:2 .•Pl /HAVEL 
AMIN  = WAVEN  • RmIN 
amax  = wavEn  * RMAX 
WRITE(6,20)  AMIN, amax 

?0  F ORMAT( /,  'ALPHA  RANGE  =',2E7.2) 

C 

c search  leg  file  for  data 

C 

REWIND  11 

rewind  12 
LESS=0 
1=0 

25  REAr(ll,END:50)  A.P.SK.L 
IF (A .GT .AMAX  ) GO  TO  50 
IF(A.L T.AMIN)  GO  TO  35 
IF(P.NE.PM)  GO  TO  35 
1 = 1+1 
ALPHA! I ) s A 
XY7  = A 

alphfr=a 
SCKRO( I )»SK 
I I ( I )=L 

IF(L.GT.LESS)  LESS=L 
READdl)  (C0EF(MM),MM  = 1,L  ) 

WRITE! 12) (C0FF( MM  ) .MM  = 1 ,L  ) 

GO  TO  25 
75  REAP! 1 1 ) 

GO  TO  25 
50  NRS  = I 

M0RrER=2. *ALPHER+3. 

REWIN1  11 
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- --  ’ M'Tkif  iT  i"  I I 


non 


REMlwn  12 


determine  Particle  pistribution  at  data  points 


no  60  I=1,NRS 
RI ( ! )=ALPHA( I )/WAvEn 
60  CONTINUE 

IF(NRS.GT.1)Go  to  6? 

OeLRI 1 ) =RMAX-RM  IN 
Go  TO  72 

62  DO  70  1=1, NRS 

IF( I .GT.I.AND. I .LT.NRS)  GO  TO  110 

IF( I . e0.NRS)DELR( I )=RMAX-(RI ( I-l)+R I ( I ) )/2, 

IF(  I .EO.l  )DELR(  n = (RI  ( I+l  )-*-RI  ( I ))/2  .-RMIN 
GO  TO  70 

no  l)ELR(  I)=(RI  (I+1)-RI  (I-l))/2.0 
70  CONTINUE 
72  WRITE(6,I12) 

n2  FORMAT!/.  'PARTin  E DISTRIBUTION  PARAMETERS') 
no  IF(NEQ)120. 132,1-10 
1.72  REAn(5.620)Pl.P?,P3,P4 
WRITE(6,630)P1,P2,P3.P4 
no  137  1=1, NRS 

1 .7  7 ENR(  n^Pl^RI  ( I )**P2«EXP(-p3*RI  ( I )»*P4) 

GO  TO  lijO 

]4C  REAn{5.620)Pl,P2,P3,P4 
WRITE(6,640)P1,P2.P3,P4 
DO  118  1=1. NRS 
IF(RI ( I ) .GT.P2)  GO  TO  l57 
ENR( I )=Pl 
GO  TO  118 

'57  ENP( 1 )=P3*RI ( I )»*P4 

118  continue 

GO  TO  150 
120  WR!TE(6.125) 

125  FORMAT!/, 'DO  NOT  KNOW  WHAT  TO  DO  WITH  NEGATIVE  NEO.  STOP.') 
STOP 

620  FoRMAT< 4FiO  .0  ) 

630  FORM*  r(/,2X,'Pl  = • ,1  PE15 .6  ,/ ,?X, 'P2  = ' ,1PE15  .6  ,/ .2X. 

1 'P3  = ',lPE15.6,/.2X,'P4  = ',lPE15.6) 

640  FORMAT! /,2X, 'PI  = ' , 1PE15. 6. /,  'P2  = ' , 1 PE  15 . 6 , / ,2 X , 

1 'P3  = '.lPE15.6./.2X,'P4  : ',iPE15.6) 

150  WRITE!6,200  ) 

200  FORMAT!/, 'DATA  POINTS '/,5X. 'ALPHA ' , 8X , ' R ' , 8 X . ' DELR ' ,13X, 'ENR' ) 
WRITE!6.220) ! I .ALPHA! I ) , R I ! I ) . DELR ! I ) ,ENR! I ) . 1=1, NRS) 

220  FORMAT! I3.F7.2.rin.2.2El5.6) 

calculate  mean  macro  x-s  values 

FN=0. 0 

no  250  1=1. NRS 

SIG!I )sSCKRO! I ) •E NR ! I ) *PI •R I ! I ) -R I ! I ) •DEL R( I ) •! .OE-08 
ENsEN-^ENR!!  )»DELR!  I ) 

L=LI  ! I ) 

READ!12)!C0EF!MM) ,MMsl,L) 

WRITE! 6, 225) I , ALPHA! I ) .SCKRO! I ) .SIC! I ) 
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o o rj  o rs 


?55  FORHaK/. 'DATA  FOR  NO  . ' » 1 3 . / . ' ALPHA  » • »F  7 . ? , 3X  . ' SCKR0»  ’ . 

1 E15.6.iX,*SlG='.E15.6) 

WRI TE (6 .230 ) ( COEF (MM) .MM=1»L) 

FORMAK/.  'COEFFICIENTS'  ,/ , ( lP6El5 .6  ) ) 

00  240  MM=1.LESS 

C0EFXS(MM  ) = C0EFXS(MM)*C0EF(MM  )*SIG(  I) 

?4P  rONTINUE 
’50  CONTINUE 
C 

r WRITE  RESULTS 
C 

WRI TF (6 ,290 )EN 

290  rORMATC/, 'MACRO  X-S  COEFF 1 C I ENTS ' , / . ' PART  I CLE  DENSITY  = *£15.6) 
WRI TE (6,3  00 ) (COErxS(MM) ,MM=1,I  ESS) 

300  FORMAT!/, {1P6E15. 6) ) 

other  output 

READ(5,302) lOUTl, I0UT2, IOUT3, I0UT4 
3(12  FoRMAT(4(1X,A4)  ) 

CHECK  HOLLERITH  VALUES 
C 

WRITE(6.303) lOUTl , 10UT2, I0UT3. I0UT4, I AO , I A1 , I A2 , 1 A3 , I A4 
303  FORMAT (/, 'OPTIONS  CHOSEN:  * , 4 ( 1 X , A4  ) , / , ' A V A I L ABLE  OPTIONS:', 

1 5(1X.A4)) 

CCCCCCCCCCCCCCCCCCCCCCCCC 

C 

C 15=1  IF  X-S  SET  FOR  ANISN  (OR  OTHER  OAK  RIDGE  CODE)  PRODUCED 
IF ( IOUTi.EQ . I AO .OR. 10UT2.EQ. 1 AO .OR. IOUT3.FO . I AO . 

\ 0R.I0UT4.EQ.  lAO)  11=1 
I iWASrl 1 

C 11=2  IF  X-S  SET  FOR  DTF  (OR  OTHER  LOS  ALAMOS  CODE)  PRODUCED 
1F(  lOllTl.EO.  lAl.OR.  I0UT2.EQ.  lAl.OR.  I0UT3.E0.  lAl. 

1 OR. I0UT4.EQ. lAl ) Il=2 
C 12=1  IF  SMOOTHING  TESTED 

IF ( IOUTi.EQ. I A2. OR. I OUT  2 . EQ . I A2 . OR . IOUT3 . EO . I A2 . 

1 0R.I0UT4.EQ. 1A2)  12=1 
C 13  = 1 IF  REFITTING  OF  COEFS  TO  F TESTED 

I F( IOUTi.EQ .1  A3 .OR. IOUT2.FQ . 1 A3 .OR. I OUT 3. EO. 1 A3 . 

1 OR . IOUT4 .EQ. IA3)  13=1 

c 14=1  If  truncation  Tested 

IF( IOUTI.EQ . I A4 .OR. IOUT2.EQ . I A4 .OR. 10UT3.E0 . I A4 . 

1 OR.IOUT4.EQ. IA4)  14=1 
C 

1F( 11 ,EO. 0)  GO  TO  288 
I F( I IWAS. EQ .1 . AND. 1 1. FO.2 ) GO  TO  270 
GO  TO  288 
220  11=1 

WRITF(6,272) 

2/2  FORMAT (/, 'CANNOT  PRODUCE  O.R.  AND  L.A.  FORMAT  X-S  ' 

1 'SETS  IN  THE  SAME  RUN.  O.R.  (E.G..  ANISN)  CHOSEN.') 

287  FORMATC  output  INDICES! *. 15, 5X  .315) 

289  FORMAT!'  OUTPUT  INDICES:  '.!iX.4l5) 

'88  IF(I  ] .EO.O)  GO  TO  305 
IFdl.F.O.DGO  TO  285 
no  280  MM=1,LESS 

rOEFXS(MM)=COEFXS(MM)/(2.*  FLOAT ( MM )♦!. ) 


nooooo  noon 


7«0  CONTINUE 

MRITE(6.289) 11, 1?. 13. 14 
GO  TO  2Bl 

WRlTE(6,2a7)ll, 12, 13. 14 
2B1  ZERO=0.0 

WRlTE(t‘>,4)  (T1TLE(  1 ).  1=1.20) 

WRI  TE(1?».291)RM1N,RMAX.PM,WAvEL,NE0.P1.P2.P3.P4,NRS,LESS.EN,M0RDER 
291  FORMAT (/, 'PARAMETERS' , / , 5E1 5 . 6 . / , 1 4 , 4E15 . 6 . / . 2 1 5 . E15 . 6 . 1 5 ) 
WR1TE(15,292)ZER0.ZER0.C0EFXS(1)  .COEFXSd) . 
t (ZERO. ZERO, ZEHO. COEFXSCMM) .MM=2. LESS) 

?9?  F0RMAT(/, 'aElZ.S  FORMAT  CROSS-SECT  1 ONS 14*» /2F12 • I . 

1 2Ei2.5,/,(3F12.1,E12.5)) 
ir( I1.EQ.1)WRITE(6,293) 

290  roRMAT(/, 'MACR  X-S  PATA  WRITTEN  ON  MAC  FILE  (l5).-  OPTION  XSOR') 
IFdl.EO.?)  WRITF(6.294) 

294  F0RMAT(/, 'MACR  X-S  DATA  WRITTEN  ON  MAC  FlLEdS).  - OPTION  XSLA') 

CALL  TESTAN(0,COEFXS.LESS) 

STOP 

end 

subroutine  TESTAN( I .AN.LESS) 

DIMENSION  X(181 ) . AN(500) . AN?(500) . AN3(500) .P(500) , 

1 P2 (500  ) .P3(101  ) .F( 181) ,F2( 181) .F3(181) ,N(181 )»  TW(181 ). 

2 FFlTdSl) 

COMMON  ALPHER.PI . IP, IP2. IP3 
COMMON/OUT/Il , 12, 13,14 

IF(LESS.LT.501)  GO  TO  10 
WR1TE(6,5)  less 

5 FORMAT(/‘, 'LESS  = '.M,'  MljST  CHANGE  DIMENSION  iN  TESTaN  AND  LEGCHK 
1'  ) 

STOP 

10  IFd.GT.i)  GO  TO  4? 

Ti=pi/iao. 

DO  30  J=l.l8l 
T=T1*( J-1) 

X( J)=COS(T) 

30  continue 

PRODUCE  leg.  polys  AND  STORE  WITH  ALL  ORDERS 
FOR  ONE  angle  ON  RECORD 

IP  = 8 
IP2  = 9 
REWIND  IP 
REWIND  IP2 
no  41  J=l,181 

P(1  ) = 1. 0 
P(2)=X( J) 

DO  40  L=3,LESS 
FLsL 

P(L)=((2.0*FL-3.0)*(X(J)*P(L-1) )-(FL-?.0) 

1 *P(L-2) )/(FL-1.0) 

RECURSION  WITH  SMALLER  ROUND-OFF  ERROR  TESTED: 

X J=X( J) 

PL1=P(L“1) 

PL2sP (L-2 ) 

P(L)=XJ*PL1-PL2tXJ*PL1-(XJ*PL1-PL2)/(Fl-1.0) 


c 

c 

c 


DID  NOT  APPRECIABLY  REDUCE  REFIT  COEFFICIENTS  PAST 

?*alpha+2*  therefore  not  due  to  round-off 
-10  continue 

WRITE(IP) (P(L).L=1.LESS) 

41  continue 
REWIND  IP 


C 

C 


PRODUCE  leg.  polys.  AND  STORE  WITH  ALL  ANGLES  FOR 
ONE  order  on  record 


IP3=10 
Rewind  ip3 

DO  135  Jsl.lSl 
!.  3?>  P ( J ) = 1 ■ 0 

WRITE(1P3)(P( J) .J=1.181) 

DO  136  J=l.l8l 
P2(J)  = P{J) 

136  P( J ) = X( J) 

WRITE(IP3)(P( J) .J=1.181) 

Do  141  LS3.LESS 
FI.  = L 

DO  140  J=l.l8l 
P3(J)sP2(J) 

P2( J)=P(J) 

P( J)s( (2.«FL-3. )*(X( J)*P2( J) )-(FL-2. ) *P3 ( J ) ) / ( FL-1 . ) 


C 

C 

C 

C 

C 

c 

c 

c 

c 

c 


SECOND  VERSION  OF  ALTERNATIVE  RECURSION: 


XJ=X( J) 

P2J=P2( J) 

P3J=P3( J) 

P( J)=XJ*P?J-P3J*XJ*P2J-(XJ*P?J-P3J)/(FL-1 .0 ) 


EXPECT  RESIDUALLY  LARGE  VALUES  OF  HIGH-ORDER  REFITTED 
coefficients  due  to  insufficient  INFORMATION  IN  FUNCTION 

140  CONTINUE 

WRITE(IP3)(P( J) ,J=1.18l) 

141  CONTINUE 
REWIND  IP3 
GO  TO  46 

4?  IFILESS.LE.LESWAS)  GO  TO  46 
LIsLESWAS*! 

DO  45  J=1.181 

READ( IP) (P(L)  .Lsl.LESwAS) 

Do  44  L=L1.  less 

fl=l 

P{L)=( {2.0*FL-3.0)»(X( J)*P(L-1) )-(FL-2.0) 

1 •P(L-?) )/(Fl-1 ‘O) 

44  CONTINUE 

WRITE! IP2)(P(L) .L=1*LESS) 

45  continue 
itEmp=ip 

IPs IP2 
IP?sI TEMP 
REWIND  IP 
rewind  IP2 
LWMsLESMaS-2 


iiii  ~ ■ *iTii  I 
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DO  144  1*1. LWM 
1 44  REAn(  IP3) 

REAPf lP3> (P2( J) .^  = 1 .181  ) 

READ(  IP3)  (P(J),  J=l.  181) 

DO  148  L=l1.LPSS 
FL-1 

DO  146  J=1.191 
P3(J)=P2(J> 

P2( J)=P( J) 

P( J)={ (2. •FL-3. )» (X( J)*P2( J))-( PL-2. )«P3( J) )/(PL-l 
146  continue 

WRITEnP3)  (P<  J)  .J=1.181) 

148  CONTINUE 
REWIND  IP3 
46  LESWAS  = less 

CALL  CALCP(AN.LESS.F) 

Ifd.GT.O)  RETURN 
C 

C EVALUATE  F USING  TRUNCATED  SERIES 


IF(  M.EO.DCALL  LEGCHK  ( AN . LESS ) 


C USE  SMOOTHED  VALUES  OF  F TO  GFT  I OWER  ORDER  FIT 


ir( I2.EQ.O)GO  TO  1 71 

CALL  SMOOTH(F.N.TW.FFIT.NFIT) 

L=1 

163  SUM*  0.0 

READ( 1P3) (P( J) . J=1 , 181) 

00  165  J2=1.NFIT 
NPT=N(J2) 

XJ=X(NPT ) 

S1NJ=S0RT{1 .-XJ»XJ) 

SUM=SUMtSINJ*FFIT ( J2)*TW( J?)*P(NPT) 

165  Continue 

ELM=L-1 

SUM=SUM*{?. *ELM+1 . ) «2 . »P1 
AN2(L )=SUM 

L*L4.1 

IF(L.LE.lESS)  go  to  163 
rewind  IP3 

WR1TB(6.i70)(AN2(L) ,L=1,LESS) 

1/0  FORMATT/. 'COEFS  FROM  SMOOTHED  F • . / . ( 1P6E15 . 6 ) ) 
call  rALCF(AN2.LESS.F2) 

IF{ I4.EQ.1)CALL  LEGCHK( AN2.1ESS) 

r 

C 

C INTEGRATE  F OVER  ALL  ANGLES.  RECALCULATE  AN(L)  VALUES 
C (F  NOT  SMOOTHED) 


1^1  IF( 13 .EQ. 0)RETURN 
L=1 

WRITE(6.62) 

6P  FORMAT!/. ‘RECALCULATEn  COEFFICFENTS 
1 '(SIMILAR  TO  1 DEG  SMOOTHING.)') 

63  SUMtO.O 

RFAp( IP3) (P( J) . Jsl .181) 

DO  65  J«1,180 


OPTION  AOIN.'./ 


on  onnon  non  ono 


C0SJ=(X(J)+X( J+i) )/2. 

S!NJ=SQRT(l .-C0SJ*C0SJ) 

rj=(r(j)+F(j*i) )/2. 

F:LM  = L-1 

PWT=(?. *ELM+1 . ) »(P( J)+P( J+1 ) )/2. 

SUM  = SUM«-SINJ*rj*Tl8PwT 

CONTINUE 

SUM=  «:UM*2.*PI 

WR1TE(6,70)L.SUM,AN(L) 

/O  f ORMATI 'L ,SUM,AN(L)  = • , I 5 . 1P2E 15 . 6 ) 
AN3  (L  ) = SUM 
iF(I.NE.O)  RETURN 
L=L4l 

IF (L.LE -LESS)  GO  TQ  63 

RFWIND  IP3 

r-ALL  CALCFI  AN3,LESS) 

TRY  TRUN-ATING  THE  RECALCULATED  SERIFS 


IF( 14.EQ.1 )CALL  LEGCHKI AN3. LESS ) 

RETURN 

END 

subroutine  LEGCHKIAN.I ESS) 

COMMON  ALPHER.PI. IP,IP2.IP3 

DIMENSION  Fl0<181 ) .LORDERI 7 ) , AN(500  ) .LTEST(6) .P(500  ) 
CMFrK  for  lower  LEGENDRE  ORDER  ACCURACY 


WR1TF(6,2)LESS 

2 FORMAT!/, 'TRUNCATE  TO  LESS  THAN 
1 'TEST  FOR’ ) 


14,  • 


TERMS.  (OPTION  TRUN)  ’ ,/, 


I 


FIND  TRUNCATION  TEST  POINTS 

DO  5 1 = 2,6. 

^ LTESTd  ) = 1 
ITEST=1 

I TEST (1 )s2. •ALPHEH+3. 

WRITE(6,8)LTEST(1) 

8 FoRm4T(I5,'  ( 2*ALPHA+2) ' ) 

'.OOK  FOR  LOCAL  MINIMUM 

L2=LESS/2 
6 DO  10  L=L2»LESS 

IF<  AN(L  ) .QT  . AN(l.-l)  .and.  AN(L>  .GT.  AN(L-2)  ) GO  TO  12 
10  CONTINUE 
GO  TO  15 

12  itest=itest*i 
LTEST( ITEST  )=L 
HRITE(6,13)L 

13  F0RMAT(I5,’  (local  MINIMUM)’) 

IF( ITEST.EO.5)  GO  TO  15 

L2  = L’-1. 

GO  TO  6 


c 

C LOOK  FOR  FIRST  COEFF 1 C I ENT . L T . AN < 1 ) 

r 

DO  20  i.=2,lESS 
IF(AN(L).LT.AN<1M  go  to  22 
?n  CONTINUE 
V?  ITESTrlTEsT+1 
LTEST(1TEST)*L 
WRI TF (6,P3)i 

format(I5.'  (Lt. first  COEF)’) 

c 

C TEST  points  IN  ORDER 

C 

1F( ITEST.EO.I ) GO  TO  27 
I loop=o 

ITESTM= ItEST-1 
TF(  UOOP.GT  .15)  STOP 
no  26  I=1,1TESTM 

IF(lTFST( I ) .LE.l TEST(  I + l)  ) GO  TO  26 
I temp=ltest  ( I ) 

LTEST(  I )aLTEST(  l + D 
LTESTd-iDsITEMP 
II OOP=ILOOP+l 
GO  TO  2A 

?6  continue 

C 

C CALCULATE  RESULTING  PHASE  FUNCTIONS 
C 

57  no  2S  1=1,ITEST 

LORDER  (I*1)  = LTEST(  I) 

28  CONTINUE 
C 

DO  125  Jsl»l8l 
1 25  FLO( J)=0.0 
C 

DO  140  L0=1. ITEST 
L1=L0R0ER<L0)+1 
L?=L0RDER(L0+1 ) 

DO  155  Jsl.iei 
READ! IP) (PIL) .L=l .LESS) 
no  150  L=L1.L2 

FLO( J)=FL0( J)+AN(L)*P(L)/( 4.*PI ) 

150  CONTINUE 

155  CONTINUE  g 

REWIND  IP  P 

C I 

WRITE(6,156)L2* ( FlO ( J ) . J= 1 • l8l ) I 

156  FORMAT!/. 'L  = ' . I 4 , / , IPEl 5 . 6 , / , ( 1 P6El 5 . 6 ) ) I 

1 lO  continue  I 

C I 

return  I 

end  I I 

subroutine  CALCF(AN,LESS.F)  I 

COMMON  ALPHER.PI , IP, IP?, IP3  j 

DIMENSION  P(500  UF{iei)  ,An(500)  J 

UO  55  J=l,181  j 

F(J)sO-0  , 

REAO(IP)(P{L) .L=1 .LESS)  | 


I 


DO  50  L=l,LESS 
r(.])  = r(  J)*AN<L)»f’<L  ) 

50  CONTINUE 

F(J)=r( J)/(4.«PI) 

55  continue 
REWIND  IP 

WRITE<6,60) (F(J), J= 1,181) 

60  F0RMAT(/,'F  FROM  I EG.  POLYS.,  0 TO  180  DEGREES' ,/ ,1PE15 .6 , 
1/.(1P6E15.6)) 

RETURN 

end 

subroutine  SMOOTH ( F .N.TW.FF I T.NFIT) 

DIMENSION  F(181) ,N(l8l),TW(l8t),FFlT(l8l) 

Tl=3. 1415926536/180. 

WRITE(6,?45) 

2^5  FORMATT/, 'SMOOTHING  PARAMETERS') 

READ! 5» 246 )THRESHF, WIDTH 
•>46  FnRMAT(2El0.0) 

WRI TE |6,247)THRESHF,WIDTH 
■’4  7 FORMAT!/,  ' THRESHF  , W 1 DTH=  ' ,1P?E15.6) 
no  248  J=2,90 
F1 =F( 1 )*THRESHF 
IF(F( J) .LT.F1  ) GO  TO  249 
^48  CONTINUE 

WRITE(6,250) 

'->50  FORMAT!/, 'NO  VALUE  OF  PHASE  FUNCT  IN  FIRST  90  DEG  LESS  THAN 
1 •F(i)*THRESHF  TO  BEGIN  SMOOTHING') 

60  TO  61 
^49  JWIDE=WIDTH 

wi dth=jwide 

NINT=180. /WIDTH 
DO  251  Isl.NINT 
JLAST  = 1 *JWlDE-fl 
IF(J,LT. jLAST)  GD  to  255 

251  CONTINUE 
WRITE!6.?52)  J 

252  format! /. 'THRESHOLD  ANGLE*' H,'  NOT  FOUND') 

GO  To  61 

255  NFIT=J*1-*>NINT-I 
WRITE!6,253) 

253  FORMAT!/, 'REP.  INDEX  WIDTH  F') 

N!1  )sl 

TW( 1 ) =T 1 /2. 

FF1T(1)=F(1) 

DO  256  K=2'J 
NIK  ) = K 
Tw(K)sTi 

256  FFIT!K)=F(K) 

WRI  TE  (6  ,254  )(  N!  K)  ,TW(  K)  ,FF1  T!  K)  ,K=1  ,J  ) 

254  F0RMAT(15,2E15.6) 

N! J + 1 )= ! JlAST»J)/2 
FFI TJ=0 .0 

JM=  Jl.  aST-1 

IF!N!  J-^1)  .EQ.  J)N(  J+l)sj4i 
IF! J.EQ.JM)G0  to  357 
no  257  J2=J^t.JM 

257  FFI  TUs^FITJ-^F  (J2) 
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^57  FFITJ=rFlTJ+F( JLAST )/2. 

P0INTS=FL0AT( JLAST-J)-0.5 
FFIT( J+l)=FriTJ/POINTS 
TW( j+i) sPOINTS*Tl 

WRI TE( 6.254 ) N( J+1) ,TW(J+1 ) .FF1T( J+D 

JMl  D = Jl.AST-JWlDE/2 

DO  260  K=J+2.NF IT 

JMII)=  JMID+JWIDE 

N(K)rjMlD 

jfirst=jlast 

JL  AST=JFirST*  JW  lOE 

FF1TK=  (F(JFIRST)+F(JlAST))/2. 

1F( JWIDE.EQ.l)  GO  TO  259 

JP=jFiRST+l 

JM-JL  AST-1 

DO  258  J=JP,JM 

FF  1TK  = FF  lTK-*-F(  J' 

258  Continue 

'^59  FFIT(K)=FF1TK/FL0AT(JWIDE) 

Tw<K)=Ti*FLOAt( JWiOe) 
WR1TE(6,?54)N(K) ,TW{K),FF1T{K) 

?60  continue 


RETURN 

ENf? 


APPENDIX  D.  MACEDIT  UTILIZATION  INSTRUCTIONS 


The  MACEDIT  code  produces  a normalized  set  of  macroscopic 
Legendre  coefficients  from  the  values  on  the  MAC  file.  All  values 
are  normalized  by  the  factor  necessary  to  change  the  first  coefficient 
to  some  Input  value.  The  first  value,  A^,  Is  the  total  scattering  croE 
section  and  can  be  used  as  a basis  for  manipulating  the  entire  set. 

As  examples  of  such  use,  an  approximate  set  of  A^  values  might  be  pro- 
duced from  an  Incomplete  LEG  file  or  a MIELEG  run  with  a small  number 
of  PARTS.  This  may  be  sufficient  detail  to  produce  relative  values 
(and  determine  the  order  necessary) , but  Insufficient  for  absolute 
magnitude.  Separate  information  of  the  magnitude  of  the  scattering 
cross  section  (as,  for  instance,  from  a MIE-2  code  run)  could  then 
produce  a complete  macroscopic  Legendre  coefficient  set  for  use  in  a 
transport  code. 

D.l  MACEDIT  Input  Description 

COEFO  (E15.0) 

MACEDIT  Input  consists  of  a single  value,  COEFO,  and  the  MAC 
data  file.  The  value,  COEFO,  is  the  magnitude  of  the  scattering  cross 
section  to  which  the  first  coefficient  Is  normalized  and  thus  forms 
the  normalizing  factor  for  all  other  coefficients. 


D.2 


MACEDIT  Output  Description  - Printout 


The  MACEDIT  printout  consists  of  all  of  the  parameters  found  on 
the  MAC  file  as  well  as  the  unnormalized  and  normalized  coefficients 
written  on  the  MAC2  file. 


D.3 


MACEDIT  Output  Description  - Files 


The  MACEDIT  code  reproduces  the  MAC-flle  coefficients  on  the 
MAC2  file,  then  writes  the  normalized  coefficients  on  the  MAC2  file 
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as  well.  Since  any  MACRO  r\in  or  MI’='LEG  run  (with  XSOR  or  XSIA 
options)  will  write  over  the  MAC  file,  the  MAC2  file  can  serve  as  a 
preserve  of  the  MAC  data.  (In  general,  any  MAC  or  MAC2  data  which 
Is  to  be  saved  Indefinitely  should  be  copied  to  some  other  file 
to  avoid  Inadvertent  loss.) 


MACEDIT  FORTRAN  Listing 


AO-A046  614  RADIATION  RESEARCH  ASSOCIATES  INC  FORT  WORTH  TEX  F/G  20/6 

STUDIES  IN  APPLICATION  OF  DISCRETE  ORDINATES  TRANSPORT  METHODS  — ETC(U) 
SEP  77  D LINDSTROM  F08606-77-C-0008 


C MArEDiT  - code  TO  NORMALIZE  MACROSCOPIC  LEGENDRE  COEFFICIENTS  ON  THE 
C MAC  FILE  TO  A KNOWN  SCATTERING  CROSS-SECTION  (FIRST  COEFFICIENT) 

C 

complex  pm 

DIMENSION  COEFXS( 500) .title (20) .C0EFX2r 500) .STUFF (8) 

zero  =0.0 

rewind  15 
REAn(5. 10)COEF0 
'0  FoRMATtElS.O  ) 

READ(15,4)(titLE( 1) .1=1.20) 

WR1TE(16,4)  (TITLEd  ),  1 = 1.20) 

4 F0rMAT(/.?0A4) 

WRITE(6 , 104){TITLE( I ) . I =1 , 20 ) . COEFO 
lO-*.  FORMAT!/, 'MACEdIT  FOR  COEFS  . OF  CASE  ; • . / . 20  A4/ 

1 • NORMALIZED  TO'  .El5.4) 

READ(15,91)RmIN,RMAX.PM.WAVEl.NEQ.P1 ,P2.P3.P4.NRS.lESS.En.MOROER 

91  FORMAT(/.lOX./.5El5.6./. I4,4El5.6,/,2I5,El5.6. 15) 
WR1TE(16,91)RM1N,RmAX.PM,WAVeL.Ne0,P1.P?,P3,P4.NrS.LESS.EN.M0P0ER 
WRITE</'.92)RMIN,RMAX 

92  FORMAT!/, IX, 'RMIN=' ,E15.6,5X, 'RMaXs ' . El5. 6) 

WRITE!*. 93)PM,WAVEL 

93  FORMaT!/,1X.  'PM='  ,2Ei5.6,5X,  •WAVEL«=’,El5.6) 

WRITE!6,9  4)NE0.P1  ,P2.P3,P4 

94  FORMAT!/, IX. •NEQ=' . I4,/, 'OISTRIB.  P ARAMETERS’ ,5X , 4E15 . 6 ) 
WRITE!6,95)NRS*LEsS,EN.M0RDER 

95  FORMAT!/. IX . 'NRS= ' , I5,5X, 'LESS=' , I5.5X, 'ENs ' ,E15.6, 

1 ' SUGG.  P-0RDER='I5) 

RE AD (15. 200)  (STuFp! J ) , I =1 , 8 ) , DUM . DUM . COeFXS ! I ) . DUM , 

1 (DUM. DUM, DUM. COEFXS! I ) , 1=2. LESS) 

20  0 forma T!/,7A4,/,A<1,/,  (4E12.5)  ) 

WRI TE( 16. 201) (STUFF (I ). 1=1,8) .ZERO. ZERO .COEFXS! 1) .COEFXS! 1 ) . 

1 (ZERO. Zero, ZERO. COEFXS! I ). 1=2. LESS) 

901  F0RMAT!/,7A4./, A4,/,2F12.1.2Et2.5,/, (3F12.1,E12.5)  ) 

FACTOR  = COEFO/COEFXSd ) 

DO  210  1=1, LESS 

C0EFX2( I )=FACTOR»COEFXS( 1 ) 

210  CONTINUE 

WRITE (6. 220) (STUFF! 1 ), 1=1,8) 

220  FORMAT!/, tX,7A4,/, IX, A4) 

WRITE (6, 240) (COEFXS! I ) ,C0EFX2! I ) , 1=1. LESS)  . 

.'40  FORM AT! /.IX, 'COEFFICIENTS' ./,7X, 'ORIGINAL' .7X, 'NORMALIZED'/ 

1 !1X,2E15.5)) 

WRITE!1.6,250)COEFO,STUFF!8)  .ZER0,7tR0,C0EFX2!l)  ,C0EFX2!1 ) . 

1 (ZERO, ZERO, ZERO, COEFX?! I ) , I =2, LESS) 

FORMAT! /, 'COEFFICIENT  NORMALIZED  TO  ' , E15 . 6 , / , A4  , / , 2F  12 . 1 , 2El2 . 5 , 

1 /. (3F12. 1,E12.5)  ) 

STOP 

END 


APPENDIX  E.  TDA  UTILIZATION  INSTRUCTIONS 


(From  CCC-180,  Time-Dependent  Multigroup  One-Dlmenslonal 
Discrete  Ordinates  Transport  Code.  RSIC,  ORNL) 


E.l  Source  and  Initial  Conditions 

Time-dependent  ANISN  offers  a choice  of  two  types  of  sources 
and  one  initial  condition  specification.  A space  and  energy  distributed 
source  with  a step  function  time  distribution  is  available.  This  source 
is  set  equal  to  zero  by  the  program  after  the  first  time  interval.  An 
analytic  first-collision  source  provides  an  accurate  representation  of 
delta  function  or  time  dependent,  point  (in  space),  isotropic  sources 
in  spheres,  infinite  plane,  mono-directional  sources  in  slab  geometry, 
and  infinite  line,  isotropic-perpendlcular-to-the-llne  sources  in 
cylinders.  All  sources  may  have  an  arbitrary  distribution  in  energy. 

The  complete,  centered  in  space  and  angle,  flux  distribution  at  T>0 
may  be  specified  alone  or  in  addition  to  one  of  the  above  sources.  If 
both  a source  and  the  initial  condition  are  specified,  the  normalization 
factor  should  be  zero  since  only  one  of  the  arrays  will  be  normalized. 

E.2  Output  Edit 

It  is  apparent  that  attempts  to  print  all  the  information 
generated  by  a time-dependent  discrete  ordinates  program  in  each  time 
interval  would  be  impractical  and  indeed  useless.  Time-dependent  ANISN 
generates  a tape  (or  disk,  etc.)  which  contains  the  scalar  flux,  the 
out-going  angular  flux  at  both  boundaries  for  each  time  interval,  and 
the  fission  neutron  density.  If  the  analytic  first-collision  source 
option  is  used  and  the  user  desires,  the  uncollided  flux  is  also  written 
on  this  tape.  In  addition,  the  user  may  specify  that  the  complete  angular 
flux  in  each  time  interval  be  written  on  a separate  tape.  Subroutine 
EDIT  is  designed  so  that  it  may  be  easily  modified  to  suit  the  demands 
of  the  user.  All  parameters  and  any  ANISN  data  arrays  which  could 


conceivably  be  used  In  an  analysis  of  these  data  sets  are  available 
for  use  In  the  subroutine. 

The  version  of  EDIT  distributed  with  the  program  will,  at  the 
users'  option,  print  the  scalar  flux,  the  uncolllded  flux  If  available, 
and  compute  activities  from  user-supplied  activity  cross  sections. 

E.3  Input  Specifications 

All  numerical  data  Is  written  In  the  FIDO  format  used  In  ANISN. 
A complete  description  of  the  format  and  convenience  options  Is  found 
In  Appendix  A.  Since  familiarity  with  ANISN  Is  assumed,  the  following 
data  description  Is  brief  except  for  those  parameters  or  arrays  which 
are  now  concerned  with  time  dependence.  As  In  ANISN,  the  quantity  In 
brackets  Is  the  array  dimension  and  the  expression  In  braces  Is  the 
condition  requiring  entry  of  an  array  or  set  of  arrays.  If  no  condi- 
tion Is  specified,  entry  of  the  array  or  set  la  required.  A T follows 
each  set  which  Is  entered. 

A.  Title  card  - format  (12A4,  18X,  16) 

A time  limit  In  minutes  may  be  entered  In  columns  67-72. 

The  case  Is  terminated  following  the  first  time  Interval 
In  which  the  time  limit  Is  exceeded. 


B.  Parameters 

15$  Integer  parameters  [36] 

1.  ID  problem  Identification  number 


ITH  - 


3.  ISCT 


4. 

5. 


ISN 

IGE 


0 - forward  solution 

1 - adjoint  solution 

maximum  order  of  Legendre  polynomial  approximation 
to  scattering  cross  sections 
angular  quadrature  order 

1 - slab 

2 - cylinder 

3 - sphere 
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6. 

IBL 

left  boundary  condition 

0 - vacuum 

1 - reflection 

2 - periodic 

3 - whlte/albedo 

7. 

IBR 

right  boundary  condition,  same  options  as  IBL 

8. 

IZM 

number  of  zones 

9. 

IM 

nvimber  of  mesh  Intervals 

10. 

lEVT 

0 

11. 

IGM 

number  of  energy  groups 

a 

CM 

rM 

IHT 

position  of  cross-section  table 

13. 

IHS 

position  of  0 in  cross-section  table 
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14. 

IHH 

length  of  cross-section  table 

15. 

MS 

length  of  cross-section  mixing  table 

16. 

MCR 

number  of  cross-section  sets  to  be  read  from 

cards 

17. 

MTP 

number  of  cross-section  sets  to  be  read  from 

library  tape 

18. 

MT 

total  number  of  cross-section  sets 

19. 

IDFM 

0 - no  effect 

1 - enter  density  factors  (21*) 

20. 

IPVT 

0 

21. 

IQM 

0 - no  effect 

1 - enter  distributed  source  to  be  used  in 

first  tine  interval,  only 

22. 

1PM 

0 - no  effect 

IM  - enter  complete  centered  angular  flux 

distribution  at  T-0,  (IM,  MM,  IGM) 

23. 

IPP 

0 

24. 

IIM 

inner  iteration  maximum  per  group  per  outer 

Iteration 

25. 

IDl 

0 - no  effect 

1 •>  print  scalar  flux 

2 > print  uncolllded  flux  If  IFG  0 

3 - both  1 and  2 
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26. 

ID2 

0 

27. 

ID3 

number  of  time  Intervals 

28. 

ID4 

number  of  activities 

29. 

ICM 

outer  iteration  maximum 

30. 

IDAU 

0 

31. 

IDAT2 

0 

32. 

IFG 

0 - no  effect 

N - enter  N source  spectra  (23*)  for  first- 

collision  source  calculation 

33. 

IFLU 

no.  of  time  steps  for  first-collision  source 

34. 

IFN 

1 - enter  flux  guess  (3*) 

35. 

IPRT 

0 - no  effect 

1 - do  not  print  cross  sections 

36. 

IXTR 

0 - no  effect 

1 - write  angular  flux  tape  for  each  time  Interval 

’ Floating  ] 

point  parameters  [14] 

1. 

EV 

0.0 

2. 

EVM 

0.0 

3. 

EPS 

accuracy  desired 

4. 

BF 

buckling  factor,  normally  1.420892 

5. 

DY 

cylinder 

6. 

DZ 

plane  depth  for  buckling  correction 

7. 

DFMl 

transverse  dimension  for  void  streaming  correction 

8. 

XNF 

normalization  factor 

9. 

PV 

0.0 

10. 

RYF 

0.5 

11. 

XLAL 

point  flux  convergence  criterion 

12. 

XLAH 

0.0 

13. 

EQL 

0.0 

14. 

XNPM 

0.0 

T 
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C.  Cross  Sections 

13$  cross-section  library  ID  numbers  [MTP]  {MTP>0} 

14*  MCR  cross-section  sets  [IHM*IGM*MCR]  {MCR>0} 

T 

D.  Distributed  source-centered  flux  distribution  {l(^«l  or 
17*  distributed  source  flM*IGM]  {l(^“l} 

18*  centered  angular  flux  distribution  [IM*MM*irai]  {IPM-IM} 

T 

E.  Flux  guess  {lFN~l} 

3*  flux  guess  [IM*I(»f]  {iFN-l} 

T 

F.  Remainder  of  data 

1*  fission  spectrum  [IGM] 

4*  radii  [IM+1] 

5*  velocities  [IGM] 

6*  angular  quadrature  weights  [MM] 

7*  angular  quadrature  cosines  [MM] 

8$  zone  numbers  [IM] 

9$  material  numbers  [IZM] 

10$  mixture  numbers  [MS] 

11$  component  numbers  [MS] 

12*  number  of  densities  [MS] 

19$  order  of  Legendre  polynomial  approximation  to  scattering 
cross  section  [IZM] 

20*  activity  cross  sections  [IGM*ID4]  {lD4>0} 

21*  density  factors  [IM]  {iDFM-l} 

22*  time  interval  boundaries  [ID3-I-1] 

23*  IFG  source  spectra;  source  emission  occurs  at  corres- 
ponding times  specified  by  24*  array;  the  first  spectrum 
is  normalized  to  XNF  (if  XNF>0)  and  the  other  spectra  are 
scaled  to  maintain  the  same  relative  values  [IGM*IFG] {IFG>0} 


24*  source  emission  time  interval  boundaries  for  first- 
collision  source  spectra  [IFLU+1]  {lFG>0} 

25*  right  boundary  albedo  [IGM]  {lBR-3} 

26*  left  boundary  albedo  [IGM]  {lBL-3} 

27$  source  spectrum  no.  by  time  step  [IFLU] 


E.4  Data  Array  Formats 

All  input  data,  with  the  exception  of  the  title  card,  are  read 
using  the  same  format. 

E.4.1  Type  1 Format 

Each  card  is  divided  up  into  six  12-dlgit  data  fields  which  are 
in  turn  divided  up  into  three  subfields.  Illustrated  in  the  following 
figure.  Only  one  data  field  is  shown. 


I ini 


The  first  subfield  is  a two-digit  Integer;  the  second  subfield  contains 
either  a $,  *,  R,  I,  T,  S,  F,  A,  C,  E,  Q,  L.  N.  M,  0,  U,  V,  Z,  +,  -,  or 
a blank.  The  third  subfield  contains  either  a fixed  or  floating  point 
number.  The  contents  of  the  first  two  subfields  will  define  the  opera- 
tion to  be  performed  on  the  third  field. 

Blank  fields  are  Ignored.  One  can  use  any  or  all  fields  on  a 
card.  For  example,  a box  of  blank  cards  sandwiched  anywhere  in  a data 
array  would  be  completely  Ignored. 

Each  data  array  is  identified  by  a two-digit  Integer  in  a first 
subfield.  There  are  both  fixed  and  floating  point  arrays;  a fixed  point 


array  is  designated  by  a $ in  the  second  subfield,  a floating  point 
array  by  an  *. 

The  second  sub field  contains  an  operator  which  specifies  the 
type  of  operation  to  be  performed  on  the  data.  The  possible  operators 
are  listed  below. 


E.4.2  Array  Operators 

$ Indicates  the  beginning  of  an  Integer  array.  The  first  sub- 
field contains  a one-  or  two-digit  number  Identifying  the  array. 

* Indicates  the  beginning  of  a floating  point  array.  The  first 
subfield  Identifies  the  array. 

R Indicates  that  the  entry  In  the  third  subfield  Is  to  be 
repeated  the  number  of  times  specified  In  the  first  subfield. 

I Indicates  linear  Interpolation  between  the  entry  In  the 
third  subfield  and  the  entry  In  the  third  subfield  of  the  next  data 
field.  The  number  In  the  first  subfield  gives  the  number  of  points 
to  be  placed  equally  spaced  In  the  specified  range. 

T Indicates  termination  of  data  reading  for  a block.  XLACS 
can  require  several  data  blocks  and  each  block  must  be  ended  with  a T. 

A block  can  contain  any  number  of  arrays.  Data  on  a card  after  a T 
will  be  Ignored. 

S Indicates  skip.  The  first  subfield  defines  the  number  of 
entries  to  be  skipped.  The  third  field  can  contain  the  first  entry 
following  the  skips.  A blank  third  subfield  would  be  Ignored. 

F Is  used  to  fill  the  remainder  of  an  array  with  the  Item 
given  In  the  third  subfield. 

A Is  used  to  address  a particular  location  In  an  array.  This 
location  Is  specified  In  the  third  subfield;  the  first  subfield  Is 
blank. 
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C Is  used  to  obtain  a count  of  the  number  of  Items  read  into 
an  array  up  to  the  point  where  the  C Is  placed.  An  Integer  ZZ  In  front 
of  the  C will  be  used  as  Identification  In  producing  a message  as 
follows : 

XX  ENTRIES  READ  IN  THE  YY  ARRAY  at  ZZC. 


E may  be  used  to  end  specifying  data  for  an  array.  This  option 
is  particularly  useful  when  it  Is  desired  to  replace  only  some  items 
In  a parlcular  array.  The  Items  In  question  are  replaced,  and  the  use 
of  an  E prevents  having  to  count  and  skip  to  the  end  of  the  array. 

Q Is  used  to  repeat  sequences  of  numbers.  The  length  of  the 
sequence  Is  defined  In  the  third  subfield.  The  number  of  times  to 
repeat  the  sequence  Is  given  In  the  first  subfield. 

L Is  similar  to  I except  that  a logarithmic  Interpolation  is 
performed  between  the  entry  points.  This  option  Is  particularly  useful 
for  defining  energy  structures  equally  spaced  in  lethargy. 

N Is  used  to  repeat  a sequence  of  numbers  In  reverse  order. 

The  length  of  the  sequence  Is  defined  In  the  third  subfield.  The 
number  of  times  to  repeat  a sequence  Is  given  In  the  first  subfield. 

M la  used  to  negate  and  repeat  an  Inverted  sequence.  The  length 
of  the  sequence  Is  given  In  the  third  subfield.  The  number  of  times  to 
repeat  a sequence  Is  given  In  the  first  subfield. 

0 Is  used  to  turn  on  (or  off)  the  card  image  edit  of  ANISN  input 
data.  As  with  the  C option,  an  Integer  In  front  of  the  0 Identifies 
the  particular  entry.  The  default  (starting)  condition  Is  not  to  edit 
the  data. 

U Is  used  to  replace  the  ANISN  Input  format  for  an  array.  The 
array  number  Is  given  in  the  first  subfield.  The  format,  written  in 
normal  FORTRAN,  Is  specified  on  the  card  immediately  following  the  card 
containing  a U.  The  parentheses  normally  capsulatlng  a format  should 
be  included. 
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V specifies  that  the  array  identified  in  the  first  subfield 
will  be  read  according  to  the  last  variable  format  read  in. 

Z is  used  to  specify  a string  of  zeros;  e.g.,  49Z  would 
place  forty-nine  zeros  into  an  array. 

+ or  - indicates  exponentiation.  The  data  in  the  third  field 
is  multiplied  by  10—  , where  N is  an  Integer  in  the  first  subfield. 
This  option  allows  one  to  specify  a number  in  up  to  nine  significant 
digits. 

Integer  data  in  the  third  subfield  must  be  right  adjusted. 
Floating  point  data  may  be  written  with  or  without  an  exponent.  If 
the  decimal  is  omitted,  it  is  assumed  to  be  immediately  to  the  left  of 
the  exponent  field.  If  there  is  no  exponent,  the  decimal  point  is 
assumed  to  be  to  the  extreme  right  of  the  nine-column  subfield. 

E.4.3  Input  Restrictions 

The  following  restrictions  must  be  observed  when  using  the 
ANISN  input  format: 

(1)  Blank  data  fields  are  Ignored. 

(2)  If  the  interpolation  option  (I)  is  used,  the  next  data 
field  may  not  be  either  blank  or  an  A entry. 

(3)  The  third  subfield  of  a data  field  containing  a $ or  an  * 
may  contain  an  integer,  N.  The  next  data  entry  is  assumed  to  be  the 
(W-1) th  member  of  the  array.  Normally  the  third  subfield  is  a blank 
and  is  Ignored. 

(4)  All  data  arrays  must  be  filled  with  the  correct  number  of 
entries.  A data  array  is  ended  by  either  starting  a new  data  array  or 
by  ending  a data  block. 
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The  transferral  of  input  data  to  input  forms  or  punched  cards 
for  a code  requiring  significant  amounts  of  input  is  always  a time 
consuming,  distasteful,  and  error-prone  process.  The  original  ANISN 
formats  were  designed  to  help  reduce  these  difficulties.  The  options 
are  convenience  features.  The  usefulness  of  the  "F"  option  which  fills 
an  array  is  obvious,  but  it  is  somewhat  harder  to  see  the  practical 
uses  for  some  of  the  more  obscure  ones  like  N,  M and  Q;  however, 
frequent  use  will  turn  up  situations  where  these  options  are  Invaluable 
For  example,  the  cosines  are  negated  and  reflected  about  90°,  a fact 
which  suggests  the  use  of  the  M option. 


There  are  justifiable  complaints  with  the  input  formats;  for 
example,  where  convenience  options  are  not  applicable,  data  can  be 
hard  to  write  because  of  the  manner  in  which  the  data  fields  are  spread 
on  the  card.  This  is  especially  true  of  integer  arrays,  where  the  data 
are  right  adjusted  in  12-column  fields.  The  ANISN  input  forms  help 
to  some  extent,  but  the  actual  keypunching  is  still  troublesome. 


The  input  format  has  been  greatly  Improved  by  Ward  Engle  of 
ORNL  who  has  designed  and  Implemented  an  all-FORTRAN  free-form,  ANISN 
input  scheme  which  has  data  items  separated  by  blanks  (as  others  do) , 
but  still  allows  all  of  the  Important  convenience  features  of  the 
earlier  formats.  The  restrictions  on  the  use  of  this  input  are  essen- 
tially that  the  user  write  the  data  in  a form  that  he  can  interpret 
within  the  context  of  the  ANISN  options.  Data  is  easily  written  and 
keypunched,  since  there  is  no  worry  about  which  type  character  falls 
in  which  columns  or  how  many  blanks  are  left  between  entries. 


The  free-form  input  can  be  Interspersed  with  the  fixed  form 
input.  To  select  free-form,  an  array  is  identified  as  either  a $$  or 
a **  array,  for  Integer  and  floating  point  arrays,  respectively. 


(1)  Any  third  subfield  (data  entry)  must  be  followed  by  one 


or  more  blanks.  This  Is  an  obvious  restriction,  otherwise  data  Inter- 
pretation would  be  Impossible. 


(2)  Only  columns  1-72  are  used. 

(3)  Numbers  with  exponents  must  not  have  Imbedded  blanks; 
e.g.,  use  l.OE+4,  not  1.0  E-W  or  1.0E+  4. 

(4)  The  old  + or  - options  (2nd  subfield)  are  not  operational. 
Significance  requirements  which  led  to  the  development  of  this  option 
can  be  had  directly. 

(5)  No  more  than  9 digits  In  a number  can  be  entered.  The 
exponent  Is  not  counted;  e.g.,  9234+09,  923400000+-1  will  work,  9234000000 
will  not  work.  Nlne-dlglt  accuracy  Is  clearly  beyond  the  significance 
available  for  single  precision  IBM  360  floating  point  operations. 

(6)  A blank  must  not  appear  between  items  which  fall  In  the 
first  and  second  subfields  with  the  old  format,  e.g.,  24R,  not  24  R. 

Note  that  the  99  restriction  on  the  number  of  repeats.  Interpolations, 
etc.,  has  been  eliminated. 

(7)  The  Z-entry  must  be  entered  as  738Z,  not  Z738.  The  old 
format  allowed  either. 

(8)  The  Q,  M,  N entries  must  be  specified  as  Q4,  not  4Q.  The 
old  format  allows  either.  An  entry  like  3Q4  accomplishes  the  same  as 
Q4  Q4  Q4.  This  Is  now  true  for  either  format. 

The  character  (')  In  column  1 of  a card  will  cause  the  contents 
of  the  card  to  be  listed  as  comments,  while  the  data  Is  read  In.  Column 
2 should  contain  the  proper  carriage  control  character;  e.g.,  blank, 0,1, 2, 
etc.  This  card  Is  Ignored  as  a data  card.  This  option  Is  also  available 
with  the  old  formats. 

Some  examples  of  the  new  format  are  given  below: 


1$$25R1J)_4_3Q3_2$$3R42_E_T 

The  first  25  entries  of  the  1$  array  are  I's  followed  by  0 and  4 and 


f 


i * 

I then  data  Input  to  the  array  ends.  The  T terminates  a data  block. 

I 
[■ 

42**J)  .0_0 . 166666 7_0 . 3333333_N2 
43**  -1.0  -0.8819171  0.3333333  M2 


This  example  puts  0.0,  0.1666667,  0.3333333,  0.333333,  0.1666667  in 
the  42*  array  and  -1.0,  -0.8819171,  -0.3333333,  0.3333333,  0.8819171 
in  the  43*  array. 


SAMPLE  TDA  INPUT 


TnA.  .5KM  R cloud,  . Ob-A . 05 ( MACRO ), P95 . S48 ( G ) 

81777  p 95  48  ? 3 0 1 ?7  0 1 3 4 4 0 96  0 9A  1 0 0 0 0 50  1 

0 ’ 0 0 I 1 .1  0 0 

16'=  2R0.0  .0001  6R0.0  0.5  .0001  3R0.0  T 

1 46« 


.0 

.0 

.17213E-03  .17213E-03 

.0 

.0 

.0  .43553E-03 

.0 

.0 

.0  .66264E-03 

.0 

.0 

.0  .77716E-03 

.0 

.0 

.0  .89232E-03 

.0 

.0 

.0  .99274E-03 

.0 

.0 

.0  .10876E-02 

.0 

.0 

.0  .11974E-02 

.0 

.0 

.0  .13072E-02 

.0 

.0 

.0  .14100E-02 

.0 

.0 

.0  .l533tE-02 

.0 

.0 

.0  .16179E-02 

.0 

.0 

.0  .17311E-02 

. 0 

.0 

.0  .18117E-02 

.0 

.0 

.0  .19017E-02 

.0 

.0 

.0  .19B21E-02 

.0 

. 0 

.0  .20535E-02 

.0 

.0 

.0  .21241E-02 

.0 

.0 

.0  .21952E-02 

. 0 

.0 

.0  .22539E-02 

.0 

.0 

.0  .23201E-02 

.0 

.0 

.0  .23698E-02 

.0 

.0 

.0  .24294E-02 

.0 

.0 

.0  .24708E-02 

.0 

.0 

.0  .25162E-02 

.0 

.0 

.0  .25520E-02 

.0 

.0 

.0  .25810E-02 

.0 

.0 

.0  .26074E-02 

.0 

.0 

.0  .26260E-02 

.0 

.0 

.0  .26381E-02 

.0 

.0 

.0  .26464E-02 

.0 

.0 

.0  .26506E-02 

.0 

.0 

.0  .264fl6E-0? 

.0 

.0 

.0  .26389E-02 

.0 

.0 

.0  .26305E-0? 

.0 

.0 

.0  .26037E-02 

.0 

.0 

.0  .25876E-02 

.0 

.0 

.0  .25558E-02 

.0 

.0 

.0  .25269E-02 

.0 

.0 

.0  .25374E-02 

• 0 

. 0 

.0  .25055E-02 

.0 

.0 

.0  .25022E-02 

. 0 

.0 

.0  .25278E-0? 

.0 

.0 

.0  .24955E-02 

.0 

.0 

,0  .24597E-02 

.0 

.0 

.0  .24186E-02 

.0 

.0 

.0  .23692E-02 

.0 

.0 

.0  .23098E-02 

.0 

.0 

.0  .23090E-02 

. n 

. 0 

.0  .23304E-02 

.0 

.t 

.0  .22702E-02 

.0 

.0 

.0  .22270E-02 

4 

0 39 


. 0 

. 0 

.0 

.21808E-02 

.0 

. 0 

.0 

.21517E-02 

.0 

. n 

.0 

.21130E-02 

. 0 

.0 

.0 

.20661E-02 

.0 

.0 

.0 

.19921E-02 

. 0 

.0 

.0 

.19208E-02 

.0 

.0 

.0 

.19502E-02 

.0 

.0 

.0 

.19713E-02 

. 0 

.0 

.0 

.18605E-02 

.0 

.0 

.0 

.17524E-02 

.0 

.0 

.0 

.16621E-02 

.0 

.0 

.0 

. 15724E-0? 

.0 

.0 

.0 

.14884E-02 

. n 

.0 

. 0 

.14103E-02 

.0 

.0 

.0 

.13141E-02 

.0 

.0 

.0 

.12306E-02 

.0 

.0 

.0 

.11379E-02 

. 0 

.0 

.0 

. 11784E-02 

.0 

. 0 

.0 

. 10839E-02 

.0 

.0 

.0 

.10018E-02 

.0 

.0 

.0 

.94927E-03 

.0 

.0 

.0 

.9253ne-03 

.0 

.0 

.0 

.9ie03E-03 

.0 

.0 

.0 

.93122E-03 

.0 

.0 

.0 

.92in0E-03 

. 0 

.0 

.0 

.90426E-03 

.0 

.0 

.0 

.84065E-03 

.0 

.0 

.0 

. 74183E-03 

.0 

.0 

.0 

.63671 E-03 

.0 

.0 

.0 

.55445E-03 

.0 

.0 

.0 

.41323E-03 

.0 

.0 

.0 

.34449E-03 

. 0 

.0 

.0 

.30543E-03 

.0 

.0 

.0 

.33228E-03 

.0 

.0 

.0 

.34548E-03 

.0 

.0 

.0 

.33528E-03 

.0 

.0 

.0 

.28019E-03 

. 0 

.0 

.0 

.l96eOE-03 

.0 

.0 

.0 

.14359E-03 

. 0 

.0 

.0 

. 74772E-04 

.0 

.0 

.0 

.40869E-04 

.0 

.0 

.0 

.16026E-04 

.0 

.0 

.0 

.74568E-05 

.0 

.0 

.0 

-.62533E-05 

T 

T 

3*^  FO.O  T 
1**  1.0 

4*'  ?6I0.0  5.0*04 
5*s  3. 0*10 

6 « e 

O.n  .001577  .003664  .005739  ;007790  .009808  .011785  .0l3''l3 
.0’558fl  .017389  .0l9l2i  .020773  .022337  .023808  .025180 

.02(S445  .027600  .028639  .029557  .030352  .031  020  -.031551 

.051962  .032233  .032369 

.032369  .032233  .031962  .031551  .031020  .030352  .029557 

.0'^8639  .027600  .026445  .025180  .023808  .022337  .020773 
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.0l<?i21  .017389  .015568  .013713  .011785  .009808  .007790 
.005719  .003664  .0015/7 

7 * •?> 


-.999476  -.998771  -.993530  -.984125  -.970592  -.952988 
-.931387  -.905879  -.876572  -.843588  -.807066  -.767159 
-.’24034  -.6/7872  -.628867  -.577225  -.523161  -.466903 
-.408686  -.348756  -.287362  -.224764  -.161222  -.097005 


032380 

.032380  .097005 
.466903  .523161 
.807066  .843588 
. 9*^4125  .993530 
8ST  FI 


.161222  .224764 
.577225  .628867 
,876572  .905879 
.998771 


.287362  .348756  .408686 
.67787?  .724034  .767159 
.931387  .952988  .970592 


9$!;  1 
19s:'!;  /O 

2791,0 

22'»  2610.0  1.666666-6  1.8453-6  2.2024-6  2.7381-6  3.4524 

4.3-553-6  5.4167-6  6.6667-6  8.0953-6  9.7025-6  1.1488-5  1.3452-5  1. 
23'- « 1.0 
24  B 0.0  0.0 
26»*  1.0 
27'}^  ITT 
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APPENDIX  F.  FILE  NAME  DIRECTORY 


5 


MIELEG  - FORTRAN  source  for  NIELEG  Code 

LM:MIELEG  - MIELEG  load  module 

MIELEGM  - Job  module  Co  compile,  load  MIELEG 

MIELEGR  - Job  module  to  run  MIELEG,  Including  data 

LEGEDIT  - FORTRAN  source  for  LEGEDIT  Code 
IJ1:LEGEDIT  - LEGEDIT  load  module 
LEGEDITM  - Job  module  to  compile,  load  LEGEDIT 
LEGEDITR  - Job  module  to  run  LEGEDIT,  Including  data 

MACRO  - FORTRAN  source  for  MACRO  Code 

LM:MACRO  - MACRO  load  module 

MACROM  - Job  moddle  to  compile,  load  MACRO 

MACROR  - Job  module  to  run  MACRO,  Including  data 

MACEDIT  - FORTRAN  source  for  MACEDIT  Code 
LM:MACEDIT  - MACEDIT  load  module 
MACEDITM  - Job  module  to  compile,  load  MACEDIT 
MACEDITR  - Job  module  to  run  MACEDIT,  Including  data 

LEG  - Microscopic  Legendre  coefficient  data  file 
BACKLEG  - LEG  backup  file 

NEWLEG  - LEG  production  file 

OTHERLEG  - LEG  overflow  file  for  other  refractive  Indices 

MAC  -*  Macroscopic  Legendre  coefficient  data  file,  used  to  produce  TDA 
data  files  by  MACRO. 

MAC2  - Data  file  produced  by  MACEDIT  from  data  found  on  MAC. 

i 


TDA  - Load  module  for  Time-Dependent  ANISN  jobs  requiring  COMMON 

storage  of  20K  or  less 

TDA30K  - Load  module  for  jobs  requiring  30K  or  less. 

TDA- IN  - TDA  job  data  file 

TDARUN  - TDA  run  module 

TDAMAIN  - TDA  MAIN  routine,  specifies  storage  In  CCWIMON  and  LIMl 
COMPFILE  - Job  module  to  compile  TDAMAIN 

EDITBO,  TDABO,  GUTSBO  - Other  modules  making  up  TDA  (FORTRAN  sources: 

EDIT,  TDASO,  and  GUTS  on  tape  LT#R104/.9002007E) 
TDAFILE  - Job  module  to  load  TDA 


LMANISN  - Load  module  for  steady-state  ANISN  (present  version  will  not 
accept  free-form  data) 

ANISNDATAX  - ANISN  job  data  file 
RUNANISN  - ANISN  run  module 


Angular  quadrature  data  sets  In  free-form  arrays  (6**  and  7**) 
Used  to  produce  TDA  data  files. 
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